'inn  Fti  e  cop*. 


August  1988 


UILU-ENG-88-2246 


COORDINATED  SCIENCE  LABORATORY 

College  of  Engineering 


AD-A19S  700 


MONTE  CARLO 
STUDIES  OF 
NONLINEAR 
ELECTRON 
TRANSPORT 
IN  III-V 

SEMICONDUCTORS 


Ki  Wook  Kim 


DTIC 

ELECTEI 


SEP  2  3  1988  , 


UNIVERSITY  OF  ILLINOIS  AT  URBANA-CHAMPAION 


Approved  for  Public  Release.  Distribution  Unlimited. 


8  8  9 


02  f 


a.  REPORT  SECURITY  CLASSIFICATION 

Unclassified 


OCCLASSlFICATION/OOWNO RASING  SChEOULE 

M/A 


A.  FIRFORMING  organization  report  NUMBER(S) 

UILU-ENG-88-2246 


a*.  NAME  OF  PERFORMING  ORGANIZATION  tb.  OFFICE  SYMBOL 

Coordinated  Science  Lab  utwueabu, 

University  of  Illinois  N/A 


Be.  AOORESS  (City.  Start  and  ZIP  Codt, 

1101  W.  Springfield  Avenue 
Urbana,  Illinois  61801 


•*.  NAME  OF  FUNOING/SPONSORING 

organization  joint  Services 

Electronics  Program 


Be.  AOORESS  City.  Slat*  and  ZIP  Cod*) 

800  N.  Quincy  Street 
'  Arlington,  VA  22217 


IonT^CA&.lS'  Iff tff/i f^/p^NONL INEAR  ELECTRON 


tb.  RESTRICTIVE  MARKINGS 

None 


3.  DISTRIBUTION/ AVAILABILITY  OF  REPORT 

Approved  for  public  release; 
distribution  unlimited 


9.  MONITORING  ORGANIZATION  REPORT  NUM8ERI3I 


7a.  NAME  OF  MONITORING  ORGANIZATION 

Office  of  Naval  Research 


7b.  AOORESS  (CUy.  Statt  and  ZIP  Codti 


800  N.  Quincy  Street 
Arlington,  VA  22217 


8b.  OFFICE  SYMBOL 
(If  aQglicabit) 


9.  PROCUREMENT  INSTRUMENT  IDENTIFICATION  NUMBER 


//  0  i  H  -  Sr  ~  Ct  V 


10.  SOURCE  OF  FUNOING  NOS 


PROGRAM 
ELEMENT  NO. 


PROJECT 

NO. 


12.  PERSONAL  AUTHORISI 


a,  TYPE  OF  REPORT 

'echnical 


KIHvHMKIimNllUKKI 


Kim,  Ki  Wook 


13b.  TIME  COVERED 

FROM  1985 


14.  OATE  Of  REPORT  <  Yr  .  Vfo..  Day, 

1988 


15.  page  COUNT 


cosati  cooes 


18.  SUBJECT  TgflMS  i Continue  on  reverie  i;  ncfiMry  and  identify  by  block  num b*r< 

Electron  transport  in  III-V  semiconductors,  GaAs/AlGaAs  materia 
system,  Nonequilibrium  situations,  Monte  Carlo  simulation 
method,  Analysis  of  transport  properties,  Semiclassical  Boltzmar 


k  ABSTRACT  Continue  on  ’■evert*  if  neceuory  and  iden nfy  by  bloc*  number}  transport  picture. 

Electron  transport  in  III-V  semiconductors,  especially  the  GaAs/AlGaAs  material  system,  is 
studied  in  various  nonequilibrium  situations.  Throughout  the  study,  a  Monte  Carlo  simulation 
method  is  used  for  the  analysis  of  transport  properties  in  the  semiclassical  Boltzmann  transport 
>icture.  The  present  work  essentially  consists  of  two  aspects.  The  first  topic  is  hot  electron 
-ransport  in  GaAs ,  focusing  on  the  electron  impact  ionization  effects.  The  dependence  of  impact 
ionization  rates  on  the  details  of  the  band  structure  is  investigated  by  using  two  (local  and 
lonlocal)  pseudopotential  methods.  The  spatial  evolution  of  the  ionization  rate  and  the  average 
ilectron  energy  are  studied  in  nonuniform  fields  characteristic  of  p+-n  junctions.  The  effects 
of  field  fluctuations  due  to  the  random  distribution  of  dopants  are  studied  as  well.  The  possi¬ 
bility  of  new  GaAs  electron-emitting  diodes  is  explored  numerically  and  compared  with  the  cor¬ 
responding  Si  devices. 

The  second  aspect  deals  with  the  effects  of  conduction  band  discontinuities  on  the  electron 
transport.  In  particular,  one-dimensional  heterostructures  are  modeled  to  study  the  nonlinear 
iransport  across  heterointerfaces.  First,  two  heterostructure  avalanche  (Continued  on  back) 


I 


CIS”  A I  3LTIGN.  A  V  A  .|_A8l  i_. 


49S  — 


jNL'Mi'ca  X.  same  AS  aft  __  -itic  ^SERS  _ 


21  ASSTBAC*  SECw'Ri'y  :-AS3iF  cation 


Unclassified 


OA  NAME  Or  3ESFOnSi  3lE  nOwCjal 


22b  ’E'.EPmcne  NUMBER 
,  Inctuae  Wa  Code, 


22c  0PC  >C£  SYMBOL 


1  FORM  1473,  33  APR 


SOiTiCN  OF  1  .AN  73  iS  OBSOlEt£. 


Unclassified 

_  a  -  : •_  a ••  s  1  a  • 


_ Unclassified _ 

StCUWITV  CLASSIFICATION  OE  THIS  PAQt 

*  19.  Abstract  (Continued  from  front) 


photodiodes  are  studied.  It  is  found  that  overheating,  enhanced  energy  relaxation,  and  carrier 
confinement  as  a  consequence  of  the  structure  in  real  space  have  a  pronounced  influence  on  the 
energy  and  momentum  distribution.  As  a  result,  the  energy  distribution  can  directly  reveal  the 
band  structure  of  the  material.  The  dependence  of  the  impact  ionization  rate  on  the  band  struc¬ 
tures  of  the  neighboring  layers  is  also  addressed. 

The  effects  of  a  nonequilibrium  phonon  distribution  on  the  electron  transport  are  studied  as  well 
The  phonon  distribution  can  be  considerably  perturbed  when  a  large  number  of  carriers  propagate 
across  an  interface,  experience  an  abrupt  energy  gain,  and  subsequently  relax  through  strong 
phonon  emission.  An  algorithm  is  developed  for  the  microscopic  analysis  of  phonon  dynamics.  It 
is  observed  that  the  hot  phonons  change  the  scattering  rate  significantly,  and  heat  the  electron 
energy  distribution. 

To  investigate  quasiballistic  electron  motion,  tunneling  hot  electron  transfer  amplifier  struc 
tures  are  studied  at  4.2  K.  The  numerical  results  demonstrate  the  existence  of  nearly  ballistic 
transport  in  the  base  and  in  the  collector  barrier,  and  confirm  that  the  experiments  can  indeed 
measure  the  energy  distribution  of  injected  ballistic  electrons.  The  device  characteristics, 
such  as  transfer  ratio  and  transit  time,  are  also  investigated  in  detail. 


MONTE  CARLO  STUDIES  OF  NONLINEAR  ELECTRON  TRANSPORT 
IN  III-V  SEMICONDUCTORS 


BY 

KI  WOOK  KIM 

B.S.,  Seoul  National  University,  1983 
M.S.,  University  of  Illinois,  1985 


THESIS 

Submitted  in  partial  fulfillment  of  the  requirements 
for  the  degree  of  Doctor  of  Philosophy  in  Electrical  Engineering 
in  the  Graduate  College  of  the 
University  of  Illinois  at  Urbana-Champaign,  1988 


Urh^na,  Illinois 


Ill 


MONTE  CARLO  STUDIES  OF  NONLINEAR  ELECTRON  TRANSPORT 

IN  m-V  SEMICONDUCTORS 


Ki  Wook  Kim,  Ph.D. 

Department  of  Electrical  and  Computer  Engineering 
University  of  Illinois  at  Urbana-Champaign,  1988 


Electron  transport  in  III-V  semiconductors,  especially  the  GaAs/AlGaAs  material  system,  is 
studied  in  various  nonequilibrium  situations.  Throughout  the  study,  a  Monte  Carlo  simulation 
method  is  used  for  the  analysis  of  transport  properties  in  the  semiclassical  Boltzmann  transport 
picture.  The  present  work  essentially  consists  of  two  aspects.  The  first  topic  is  hot  electron 
transport  in  GaAs,  focusing  on  the  electron  impact  ionization  effects.  The  dependence  of  impact 
ionization  rates  on  the  details  of  the  band  structure  is  investigated  by  using  two  (local  and  non¬ 
local)  pseudopotential  methods.  The  spatial  evolution  of  the  ionization  rate  and  the  average 
electron  energy  are  studied  in  nonuniform  fields  characteristic  of  p+-n  junctions.  The  effects  of 
field  fluctuations  due  to  the  random  distribution  of  dopants  are  studied  as  well.  The  possibility 
of  new  GaAs  electron-emitting  diodes  is  explored  numerically  and  compared  with  the 
corresponding  Si  devices. 

The  second  aspect  deals  with  the  effects  of  conduction  band  discontinuities  on  the  electron 
transport.  In  particular,  one-dimensional  heterostructures  are  modeled  to  study  the  nonlinear 
transport  across  heterointerfaces.  First,  two  heterostructure  avalanche  photodiodes  are  studied. 
It  is  found  that  overheating,  enhanced  energy  relaxation,  and  carrier  confinement  as  a  conse¬ 
quence  of  the  structure  in  real  space  have  a  pronouned  influence  on  the  energy  and  momentum 
distribution.  As  a  result,  the  energy  distribution  can  directly  reveal  the  band  structure  of  the 
material.  The  dependence  of  the  impact  ionization  rate  on  the  band  structures  of  the  neighbor- 


•J 


i 


ing  layers  is  also  addressed. 


The  effects  of  a  nonequilibrium  phonon  distribution  on  the  electron  transport  are  studied  as 
well.  The  phonon  distribution  can  be  considerably  perturbed  when  a  large  number  of  carriers 
propagate  across  an  interface,  experience  an  abrupt  energy  gain,  and  subsequently  relax  through 
strong  phonon  emission.  An  algorithm  is  developed  for  the  microscopic  analysis  of  phonon 
dynamics.  It  is  observed  that  the  hot  phonons  change  the  scattering  rate  significantly,  and  heat 
the  electron  energy  distribution. 

To  investigate  quasiballistic  electron  motion,  tunneling  hot  electron  transfer  amplifier  struc¬ 
tures  are  studied  at  4.2  K.  The  numerical  results  demonstrate  the  existence  of  nearly  ballistic 
transport  in  the  base  and  in  the  collector  barrier,  and  confirm  that  the  experiments  can  indeed 
measure  the  energy  distribution  of  injected  ballistic  electrons.  The  device  characteristics,  such  as 
transfer  ratio  and  transit  time,  are  also  investigated  in  detail. 
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CHAPTER  1 
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INTRODUCTION 
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Recent  advances  in  semiconductor  crystal  growth  techniques  such  as  molecular  beam  epi¬ 
taxy1,2  (MBE)  and  metalorganic  chemical  vapor  deposition3  (MOCVD)  have  made  possible  the 
fabrication  of  ultrasmall  devices,  and  have  opened  new  possibilities  for  "band  gap  engineering"  of 
semiconductor  heterostructures  with  spatially  modulated  gaps.4  By  using  the  variability  of  the 
boundary  conditions  which  can  be  imposed  on  the  electronic  wavefunction,  electrical  and  optical 
responses  in  these  structures  can  be  tailored  virtually  at  will.  Most  of  the  effort  in  this  area  has 
been  devoted  to  III-V  compound  semiconductor  materials,  especially  the  GaAs/AlxGa^As 
material  system.  Along  with  the  possibility  of  forming  latticed-matched  artificial  structures  with 
A^Ga^As,  GaAs  itself  is  known  to  have  electronic  properties  which  are  in  many  ways  superior 
to  Si.5-8  Compared  to  Si,  GaAs  shows  a  significantly  higher  electron  mobility  in  the  low  field 
regime  due  to  the  smaller  electron  mass  and  a  considerably  enhanced  drift  velocity  at  high  fields 
in  short  channel  devices  resulting  from  the  overshoot  phenomenon.  Because  of  the  relatively 
weak  momentum  relaxation  processes  at  low  electron  energies  through  polar  optical  phonon 
scattering  (in  contrast  to  deformation  potential  scattering  with  large  momentum  change  in  Si), 
electrons  in  GaAs  can  initially  gain  considerable  energy  from  the  field  without  significant  relaxa¬ 
tion  in  momentum  distribution  and  have  a  velocity  overshoot  persisting  over  a  considerable  time 
and  distance.  These  properties  and  the  direct  band  gap  make  GaAs  a  suitable  material  for  high 
speed  electronic  and  optoelectronic  devices. 

Some  of  the  devices  which  rely  on  these  properties  include  superlattice  structures,9 
modulation-doped  structures,10  quantum-well  heterostructure  lasers,11  heterostructure  avalanche 
photodiodes,12,13  real-space  transfer  switching,14  and  ballistic  transfer  structures.15,16  In  many  of 
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these  devices,  electrons  often  experience  a  reduced  dimensionality  in  the  motion  (size  quantiza¬ 
tion  effects)  and  the  displacement  in  the  momentum  distribution  due  to  the  differences  in  materi¬ 
als  between  adjacent  layers.  A  wide  range  of  novel  effects  such  as  impact  ionization  enhance¬ 
ment,12,13  real-space  transfer,14  and  hot  phonon  effects17,18  has  been  observed  in  addition  to  stan¬ 
dard  textbook  phenomena  such  as  the  particle  in  a  box  problem  (e.g.,  red  light  out  of  GaAs 
quantum  wells19),  the  Aharonov-Bohm  effect,20  and  various  tunneling  effects.21  These  effects  are 
based  on  electron  transport  far  from  equilibrium,  and  are  essential  in  the  operation  of  the  dev¬ 
ices  mentioned  above. 

Since  experimental  studies  of  the  nonlinear  transport  are  often  very  difficult  due  to  its 
microscopic  nature,  numerical  studies  with  a  microscopic  point  of  view  are  necessary  for  better 
understanding.  Conventional  device  modeling  techniques  based  on  the  drift-diffusion  equations 
are  inadequate  to  correctly  treat  the  nonlinear  transport  problems  since  these  macroscopic 
models  cannot  accurately  account  for  the  electron  motion  (or  electron  distribution)  in  momen¬ 
tum  space  which  is  essential  in  understanding  the  nonlinear  effects.  Instead,  a  technique  which 
has  shown  great  success  for  the  study  of  nonlinear  transport  is  the  Monte  Carlo  simulation 
method  22  By  stochastically  tracing  electron  movements  in  both  real  space  and  momentum 
space,  the  Monte  Carlo  method  can  generate  the  exact  numerical  solution  of  the  Boltzmann 
transport  equation  without  solving  the  complicated  integro-differential  equation  directly.23  The 
fact  that  the  simulations  are  based  upon  the  simple  first  principles  gives  the  Monte  Carlo  tech¬ 
nique  great  flexibility  in  the  inclusion  of  complicated  details  such  as  a  realistic  band  structure, 
many  scattering  mechanisms,  complicated  boundary  conditions,  and  time-  and  space-dependent 
parameters  The  many-particle  scheme  proposed  by  Lebwoh)  and  Price24  is  the  most  suitable 
when  the  collective  effects  are  appreciable  and  the  transient  information  is  of  prime  concern 
The  self-consistent  electric  field  can  naturally  be  included  in  the  simulation  by  solving  the  Pois 
son  equation  simultaneously. 
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It  is  the  purpose  of  this  thesis  to  study  the  electron  transport  in  various  nonequilibrium 
situations  in  III- V  semiconductors,  especially  the  GaAs/AlGaAs  material  system.  Throughout 
this  thesis,  a  Monte  Carlo  simulation  method  has  been  used  for  the  first-principle  analysis  of 
transport  properties  under  the  semiclassical  Boltzmann  transport  picture.  In  Chapter  2,  various 
models  to  implement  complicated  boundary  conditions  into  a  Monte  Carlo  method  will  be 
described.  The  treatments  of  coupled  plasmon/phonon  scattering,  the  Pauli  exclusion  principle, 
tunneling,  and  band  bending  will  be  discussed  in  detail.  A  quantum  mechanical  correction  (colli¬ 
sion  broadening)  to  the  Boltzmann  transport  equation  and  its  influence  on  the  electron  transport 
will  also  be  analyzed. 

Chapter  3  includes  the  study  of  hot  electron  transport  in  GaAs.  The  dependence  of  impact 
ionization  rate  on  the  details  of  the  band  structure  will  be  addressed.  The  effects  of  nonuniform 
electric  field  in  p+-n  junctions  will  be  investigated  along  with  the  fluctuation  in  the  dopant  dis¬ 
tribution.  The  results  will  be  compared  with  the  hot  electron  characteristics  of  Si  to  explore  the 
possibility  of  new  GaAs  electron-emitting  devices. 

Chapters  4  through  6  are  devoted  to  the  modeling  of  GaAs/AlGaAs  one-dimensional 
heterostructures.  The  effects  of  the  conduction  band  discontinuity  on  the  electron  transport  will 
be  discussed  along  with  the  models  used  in  the  analysis.  In  Chapter  4,  the  nonlinear  transport 
across  heterointerfaces  will  be  studied  in  two  heterostructure  avalanche  photodiodes.1213  The 
evolution  of  the  electron  energy  and  momentum  distribution  will  be  investigated  both  in  time 
and  space  The  enhancement  of  impact  ionization  in  these  structures  will  also  be  studied  in  con¬ 
junction  with  the  dependence  on  the  band  structures.  Chapter  5  will  address  hot  phonon  effects 
which  can  occur  when  a  large  current  flows  into  a  well  layer.  An  algorithm  to  incorporate  pho¬ 
non  dynamics  will  be  described 
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To  investigate  the  quasiballistic  nature  of  electron  motion  and  to  examine  the  validity  of 
hot  electron  spectroscopy,  the  tunneling  hot  electron  transfer  amplifier  structures  proposed  by 
Heiblum  et  al  ie  will  be  studied  in  Chapter  6.  Special  attention  will  be  paid  to  the  transport  in 
the  collector  barrier,  a  region  where  the  Heiblum’s  device  differs  from  the  similar  planar  doped 
barrier  transistor.15  The  switching  characteristics  will  also  be  discussed. 
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CHAPTER  2 

■ 

MONTE  CARLO  SIMULATION  METHOD 


fi 


2.1  Introduction 

The  Monte  Carlo  method  consists  of  simulating  the  motion  of  one  or  more  electrons  inside 
of  a  crystal,  subject  to  the  action  of  the  forces  (e.g.,  an  electric  field)  and  of  given  scattering 
mechanisms.  The  duration  of  carrier  free  flight  time  and  the  scattering,  events  involved  in  the 
simulation  are  selected  stochastically,  based  on  the  given  probabilities  describing  the  microscopic 
processes.  In  real  devices,  especially  heterostructure  devices,  electrons  also  experience  boundary 
conditions  such  as  reflections  at  the  interface.  The  accuracy  of  the  calculation  is  strongly  depen¬ 
dent  on  the  physical  model  of  these  microscopic  processes. 

This  chapter  is  devoted  to  a  description  of  models  used  in  the  simulation.  The  basic  Monte 
Carlo  scheme  and  the  treatment  of  conventional  scattering  processes  (such  as  electron-phonon 
interactions  and  ionized  impurity  scatterings)  have  been  well  documented  in  the  literature,25-27 
and  will  not  be  repeated  here.  Particular  attention  will  be  paid  to  more  complex  situations 
which  occur  in  heavily-doped  or  novel  structures.  In  Section  2.2,  the  treatment  of  collective 
effects,  which  are  important  for  a  degenerate  electron  gas,  will  be  discussed.  Emphasis  will  be  on 
two  topics,  coupled  plasmon/phonon  scattering  and  the  Pauli  exclusion  principle.  A  widely  used 
treatment  to  implement  tunneling,  the  transfer  matrix  method,  will  be  discussed  in  Section  2.3 
along  with  techniques  to  obtain  the  effective  mass  at  energies  within  the  band  gap.  In  Section 
2.4,  the  band  bending  due  to  the  formation  of  accumulation  and  depletion  regions  at  the  inter¬ 
faces  will  be  discussed.  The  last  section  of  this  chapter  will  deal  with  a  quantum  mechanical 
correction  to  the  semiclassical  Boltzmann  transport  equation.  The  spectral  density  technique 
will  be  described  for  the  inclusion  of  collision  broadening  into  a  Monte  Carlo  simulation. 
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2.2  Collective  Effects 

As  the  electron  concentration  increases,  interactions  between  electrons  become  nonnegligi- 
ble.  Also,  the  effect  of  band  filling  (or  the  Pauli  exclusion  principle)  which  has  been  neglected  in 
the  calculation  of  scattering  rates  cannot  be  disregarded.  Normally,  Coulombic  interactions 
between  electrons  are  categorized  into  two  scattering  mechanisms:  the  interaction  between  two 
independent  electrons  with  short  wavelength  (electron-electron  scattering)  and  the  interaction 
between  an  electron  and  a  collective  excitation  of  an  electron  gas  at  long  wavelength  (plasmon 
scattering).  A  cutoff  wavevector  qc  between  the  individual  and  the  collective  excitation  of  an 
electron  gas  is  not  well  defined.  Roughly,  one  half  of  the  Thomas-Fermi  wavevector28  has  been 
used  as  qc  where  the  collective  excitation  starts  to  damp  strongly  29  As  the  plasmon  frequency 
approaches  to  the  longitudinal  optical  phonon  energy,  the  coupling  between  these  modes  and  the 
Frohlich  polaron  (polar  optical  phonon  scattering)  becomes  significant.  In  this  section,  a  theoret- 
|  ical  model  of  these  coupled  modes  will  be  discussed.  Inclusion  of  the  Pauli  exclusion  principle 

will  also  be  explored. 


2.2.1  Coupled  plasmon/ phonon  scattering 

The  scattering  rate  of  the  coupled  plasmon/phonon  modes  can  be  calculated  by30 


(2.1) 


and 


S(q,w)  =  — Im[ 


ex(q,w) 


(2.2) 


where  ec  represents  the  dielectric  constant  in  free  space  Since  T  =  0  K  is  assumed  in  this  for¬ 
malism,  only  the  emission  events  are  considered.  Thus,  fp  is  either  0  or  I  and  is  0  m  the  cal 
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culation.  The  total  dielectric  function  calculated  in  the  random  phase  approximation  can  be 
,31 


written 

eT(q,w)  =  eL(q,o>)+€e(q,w)-l, 

where  the  dielectric  function  of  the  lattice  takes  the  form 
eL(q,W  =  l-«oo(— 2 - 2  ' 

or— a>ro 

and  the  electronic  dielectric  function  is  denoted 

£e(q,w)  =  [e0o+€reL(q>w)-1)+«  [^q.")!- 
Then,  the  total  dielectric  function  can  be  simplified  as 

6r(q,w)  -  [eoo(^=^)-Nr'L(q^)-l]+1'[fiL(q^)]- 
ur — 

The  real  and  imaginary  parts  of  the  Lindhard  function  are  written 


^  - 1+D  'OT>  & 10 1 


2_ 
1_ 

q'  2 
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I] 


^q'  2*' 


<i'L(q,w)  = 


(2.3) 


(2.4) 


(2.5) 


(2.6) 


(2.7) 


D(lSr,'K?f ' q<2kF'  0<lr“<-^2-2qk^) 

^)^Mf'^)Jl'0<_^r(q2_2qfcF)<lr"<^r(q2+2qkF) 


€0q4  ^4^q 

*_2 

0,  &o»  — -^-(q2+2qkp) 

2m 

tr2 

0,  q>2kp  tTu;< — -r(q2-2qkF) 
’  2m 


(2.8) 
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where  q'  =q/kF  and  1/  =m*S‘u'/S2kF.  The  parameter  D  is  the  number  of  equivalent  valleys  (i.e., 
1  for  the  T  valley),  and  kF  is  the  Fermi  wavevector.  When  €'L=0,  Im[l/ex(q>w)]  has  resonance 
frequencies  known  as  the  dispersion  relation  and  is  dominated  by  these  resonance  frequencies  in 
the  integration  over  the  frequency  u.  By  application  of  the  Dirac  formula32 

- 7-  .  =  P[ — i-j-]  -  )  (2.9) 

X-x'+lT]  lx~x'  v  ’ 


and  the  use  of  following  relationship 

for  the  zeros  x0  of  f(x),  the  scattering  rate  in  the  region  where  €^=0  can  be  obtained  by 

1  _  yi  r  d3q  2e2?r _ 1 _ 

’■pi-popO^)  uj  (2?r)3  f0q2  I  ^T(q^)  . 

1  dcj 


(2.10) 


(2.11) 


where  u0'  s  are  the  solutions  of  the  equation  eT(q,w)=0.  For  the  region  where  e'L  is  nonzero, 
the  resonance  frequencies  are  broadened  out  and  should  be  integrated  over  the  possible  energy- 
momentum  phase  space. 

The  resulting  dispersion  relation  gives  the  following  four  coupled  modes33  as  shown  in  the 
Fig.  2.1:  two  phonon-like  modes  at  small  and  large  wavevectors,  a  plasmon-like  mode  at  small 
wavevectors,  and  the  Landau  damping  mode  which  results  in  the  excitation  of  electrons  from  the 
Fermi  sea  through  the  damping  of  collective  motion  Depending  on  the  plasmon  frequency 
Wp={ne2/e0e00m*)1/2,  the  w+  and  u>_  branches  in  the  figure  respectively  become  plasmon-like  and 
phonon-like  for  wp>Wlo>  and  vice  versa  for  wp°<wl0-  The  Landau  damping  mode  corresponds 
to  the  region  where  e'L  is  nonzero.  In  a  different  approach,  the  Landau  damping  mode  can  be 
treated  as  individual  pair  electron-electron  scatterings  as  briefly  discussed  previously  29  To  dis¬ 
tinguish  the  collective  plasmon  interaction  and  the  individual  pair  interaction,  a  cutoff 


(meV) 


Wavevector  (q  /  k  F) 


LP- 3010 


Figure  2.1  A  schematic  plot  of  the  dispersion  relation  of  various  modes  contributing  to  the 
coupled  plasmon/phonon  scattering  processes.  A  degenerate  electron  density  of 
lX1018cm-3  and  Ej*  0.13  eV  have  been  assumed. 
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wavevector  has  been  introduced  based  upon  the  singularities  of  the  dielectric  function  in  such 
treatments.  The  approach  which  uses  the  total  dielectric  function  is  more  rigorous. 

The  calculated  scattering  rates  by  these  coupled  modes  at  n=lX1018cm-3  are  shown  in  Fig. 
2.2  as  a  function  of  energy.  For  plasmon-  and  phonon-like  modes,  the  scattering  rates  have  been 
obtained  in  the  limit  of  either  short  or  long  wavevectors.  By  making  the  short  wavevector 
approximation,  the  Lindhard  function  becomes 


er'L(q— 0,w)  =  1— €oo— ^p(q)  =  uf  (l+j— fr-q2), 


5  *2  o2 

m  wpu 


(2.12) 


which  gives  the  following  dispersion  relation: 

w!(q)  = 


2/  x  _  ^(qlWol  ±  {[u;P2(q)-F(jLo]2~~4a;p^)(*>ro^1/2 


(2.13) 


In  actual  calculations,  w±(q=0)  has  been  used.  The  phonon-like  mode  at  long  wavevectors  has 
been  treated  as  uncoupled  polar  optical  phonon  scattering  (q— K>o  approximation).  As  can  be 
seen  in  the  figure,  the  scattering  rate  is  zero  below  the  Fermi  energy  and  sharply  increases  after¬ 
wards  until  it  reaches  the  maximum.  The  Landau  damping  is  the  most  dominant  mechanism  at 
low  energies,  and  the  plasmon-like  mode  (o>+)  becomes  important  as  the  energy  increases.  The 
contribution  by  the  phonon-like  modes  (wlo  and  u;-)  is  negligible.  Figure  2.3  shows  the  coupled 
plasmon/phonon  scattering  rate  as  a  function  of  electron  density  for  an  electron  of  E=0.3  eV.  It 
is  observed  that  as  the  density  increases,  the  scattering  rate  increases  up  to  a  point  and  then 
decreases.  This  decrease  is  due  to  the  shift  in  the  Fermi  energy  toward  higher  energies  with  the 
increase  in  density,  thereby  shifting  the  onset  of  scattering.  Although  the  maximum  scattering 
rate  increases  with  density,  the  scattering  rate  with  large  density  can  be  still  in  the  increasing 
stage  and  relatively  small  at  E=0.3  eV  as  shown  in  Fig  2.4. 
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2.2.2  Pauli  exclusion  principle 

When  the  electron  gas  is  degenerate,  the  assumption  that  the  occupation  probability  f(lT)  is 
much  less  than  one  is  no  longer  valid.  The  (l-f(TT)]  factor  neglected  in  the  standard  Monte  Carlo 
simulation  is  expected  to  be  of  considerable  importance.  Inclusion  of  the  fermionic  nature  of 
electrons  (the  Pauli  exclusion  principle),  however,  requires  a  knowledge  of  the  distribution  func¬ 
tion.  Attempts  have  been  made  to  incorporate  degeneracy  into  a  single  particle  Monte  Carlo 
method  with  only  limited  success.34  The  ensemble  Monte  Carlo  method,  for  which  the  distribu¬ 
tion  can  be  obtained  naturally  at  all  times,  is  the  appropriate  choice.35 

The  procedure  generally  used  to  include  the  Pauli  exclusion  principle  is  as  follows.  The 
occupation  probability  f(lT)  is  updated  during  the  simulation  by  dividing  the  number  of  electrons 
in  each  mesh  cell  in  momentum  space  by  the  number  of  total  states  available.  Electron  drift  is 
treated  as  in  the  nondegenerate  case.  When  a  scattering  occurs,  the  final  state  is  selected  first 
without  considering  the  exclusion  principle.  Based  upon  f(lc),  the  rejection  method  is  then 
applied  to  test  if  the  chosen  final  state  is  unoccupied.  If  the  transition  is  rejected,  the  scattering 
event  is  treated  as  a  self-scattering  and  the  state  of  electron  is  not  changed. 

The  size  of  the  mesh  cells  in  momentum  space  should  be  chosen  carefully,  depending  on  the 
number  of  electrons  in  the  ensemble.  While  it  is  preferable  to  keep  the  size  of  the  mesh  small  for 
better  resolution  of  the  distribution  function,  it  must  be  large  enough  so  as  not  to  introduce  a 
significant  statistical  error  when  the  maximum  occupancy  of  a  cell  by  the  ensemble  electrons  is 
rounded  off  to  an  integer.  It  is  important  to  ensure  that  calculations  are  not  dependent  on  the 
number  of  electrons  in  the  ensemble  and  the  size  of  the  mesh  cells. 

The  method  discussed  above,  however,  requires  a  significant  amount  of  memory  in  space 
dependent  simulations.  Instead  an  approximate,  but  much  simpler,  method  can  be  used  in  dev¬ 
ice  simulations  since  the  degeneracy  effects  are  important  only  in  heavily  doped  regions  In  these 
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regions,  large  charge  density  results  in  low  electric  fields  and  low  electron  temperature.  Strong 
electron-electron  interactions  tend  to  restore  the  Fermi-Dirac  distribution.36  Therefore,  the  occu¬ 
pation  probability  function  can  be  validly  approximated  as 
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+  exp[ 


%-Erftt). 
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(2.14) 


where  the  Fermi  energy  Ep(?,t)  is  obtained  from  the  local  electron  density  and  the  electron  tem¬ 
perature  T^t)  is  from  the  local  average  energy. 


2.3  Tunneling 

As  the  control  of  epilayer  growth  becomes  more  accurate,  quantum  mechanical  effects  such 
as  tunneling  should  be  included  at  the  interfaces  along  with  classical  transfer.  In  general,  these 
effects  can  be  incorporated  into  the  Monte  Carlo  simulation  from  an  independent-particle  point 
of  view.37,38  Inclusion  of  many  particle  effects  such  as  interference39  is  inherently  difficult  due  to 
the  stochastic  nature  of  the  Monte  Carlo  method.  Although  the  effective  mass  approximation  is 
no  longer  valid  at  abrupt  interfaces,  it  is  commonly  used  due  to  the  lack  of  a  simple  alterna¬ 
tive.'40  More  precise  results  can  be  obtained  by  using  an  empirical  Bloch  wavefunction 
approach.41  In  the  calculation  of  tunneling  probability  two  methods,  the  WKB  approximation 
and  the  transfer  matrix  method,42,43  are  widely  used  under  the  effective  mass  approximation.  In 
this  section,  the  transfer  matrix  method,  which  is  the  more  accurate  of  the  two,  is  described  in 
detail. 

Figure  2.5  shows  the  schematic  of  an  arbitrary  one-dimensional  potential  barrier.  Instead 
of  using  a  continuously  varying  potential  energy,  the  barrier  is  approximated  as  a  summation  of 
constant  potential  energy  segments.  It  is  assumed  that  energy  and  parallel  momentum  are  con¬ 
served  throughout  the  calculation.  The  wavefunction  0;  in  the  jth  region  is  given  by 
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V’j(x)  =  Ajexp(j'kjx)  4-  Bjexp(— »kjx) 


(2-15) 


where  kj  is  obtained  from  the  energy-momentum  conservation  assumption.  From  the  continuity 
of  0j(x)  and  (1/mj*)  (dVj/dx)  at  each  boundary,44  determining  An  and  Bn  can  be  reduced  to  the 
multiplication  of  the  following  n  (2  X  2)  matrices: 


A0' 
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An 
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j-1 

Bn. 

Mu  M12' 
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,M21  M/J2 

Bn, 

(2.16) 


By  setting  Bn=0,  the  tunneling  probability  P(k0)  can  be  written  as 


p/i  A  (  ^n/mn  ,  |An  I  /  ^n/^n  ^  i,  /  i_2 


(2.17) 


If  nonparabolicity  is  included,  the  effective  mass  is  a  function  of  kinetic  energy  (thus,  position  in 
real  space)  even  in  the  homogeneous  material  system.  The  effective  mass  in  the  forbidden  gap 
(i.e.,  in  the  barrier  region)  can  be  defined  by  using  the  extension  of  a Tc-ff  method.45,48 

The  model  discussed  above  describes  the  tunneling  process  in  which  electrons  stay  in  the 
same  valley  in  momentum  space,  and  does  not  include  the  elastic  intervalley  transfer  between 
band  minima.  Negative  differential  resistance  has  been  reported  in  a  single  barrier  structure 
through  the  X  minima  in  the  GaAs-ALAs  system47,48  and  the  valence  band  in  the  CdTe-HgCdTe 
system.49  The  tight  binding  method  can  be  used  to  estimate  these  effects  of  the  mixing  between 
valleys.41,50 

Although  generally  used,  the  assumption  of  parallel  momentum  conservation  does  not  have 
a  rigorous  theoretical  justification  other  than  the  hand-waving  argument  of  lattice  homogeneity 
in  the  direction  parallel  to  the  interfaces.  Due  to  the  position  dependent  effective  mass,  the  total 
energy  (kinetic  plus  potential)  in  the  direction  of  propagation  is  not  conserved  even  when 
scattering  is  not  present.  This  seems  to  suggest  that  tunneling  should  be  treated  three- 
dimensionally  in  principle.  Further  study  is  needed  to  clarify  this  point. 
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2.4  Band  Bending 

In  the  tunnel  structures,  tunneling  currents  are  extremely  sensitive  to  the  actual  shape  of 
the  barrier  through  which  charge  carriers  tunnel,  and  to  the  voltage  drop  in  the  electrodes  as 
well  as  in  the  barrier.  Thus,  it  is  essential  to  perform  self-consistent  calculations  of  the  energy- 
band  diagram  which  includes  the  effects  of  band  bending  and  the  buildup  of  space  charge  in  the 
quantum  wells.  In  some  cases,  energy-band  proxies  even  indicate  that  the  dominant  carrier 
transport  mechanisms  are  quite  different  from  those  expected  when  these  effects  are  neglected. 
An  algorithm  to  obtain  realistic  energy-band  profiles  is  discussed  in  this  section,  primarily  in 
single-  and  double-barrier  structures  5,-53 

Figure  2.6  illustrates  a  schematic  of  the  potential  profile  of  a  single-barrier  structure  with 
an  undoped  AlGaAs  layer  sandwiched  between  heavily  doped  GaAs  contact  layers.  For  the  cal¬ 
culation  of  band  profile  in  this  structure,  the  following  boundary  conditions  are  generally  used: 
(i)  Far  away  from  the  barrier  layer,  the  semiconductor  properties  are  bulk-like,  (ii)  There  is  no 
interface  charge,  (iii)  Under  no  bias,  the  Fermi  level  remains  constant  throughout  the  entire 
structure  Under  an  external  bias  Va,  the  quasi-Fermi  level  in  each  contact  layer  is  constant  up 
to  the  interface  with  the  barrier.  These  quasi-Fermi  levels  must  differ  by  eVa  With  these  boun¬ 
dary  conditions,  the  Poisson  equation  is  solved  self-consistently.  Starting  from  an  initial  charge 
distribution,  the  potential  profile  and  the  charge  distribution  are  iterated  self-consistently.  Since 
the  quasi-Fermi  level  is  always  known  (conditions  (i)  and  (iii)  above),  the  charge  distribution 
which  corresponds  to  a  given  potential  profile  can  be  readily  obtained.  In  double-barrier 
resonant  tunneling  structures,  the  Schrodinger  equation  is  solved  simultaneously  with  the  Pois¬ 
son  equation  to  include  the  charge,  buildup  in  the  quantum-well  layer  The  quasi-Fermi  level  of 

the  contact  layer,  from  which  most  of  the  carriers  in  the  well  layer  come,  is  generally  used  in  (he 
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Figure  2.6 


A  schematic  of  the  potential  profile  ol  a  single-barrier  structure  with  an  undoped 
AlGaAs  layer  sandwiched  between  heavily  doped  GaAs  layers.  The  -=o!icj  line 
includes  the  band  bending,  and  the  dashed  line  does  not. 
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2.5  Collision  Broadening 

The  development  of  methods  to  calculate  the  electron  transport  properties  of  materials  far 
away  from  equilibrium  conditions  is  an  unsolved  problem  in  physics  today.  An  understanding  of 
the  behavior  of  electrons  in  high  electric  fields  is  becoming  increasingly  important  in  semiconduc¬ 
tor  structures  which  approach  angstrom  sizes.  Such  effects  as  collision  broadening,  finite  colli¬ 
sion  lifetimes,  the  intra-collisional  field  effect,  and  extremely  high  electron  scattering  rates  must 
be  considered  in  studies  of  these  nonequilibrium  problems.  For  systems  which  are  not  highly 
perturbed  from  equilibrium  and  where  such  effects  are  of  only  limited  importance,  semiclassical 
calculations  based  on  the  Boltzmann  transport  equation  prove  to  be  very  useful  in  predicting 
physical  properties.  As  mentioned  previously,  Monte  Carlo  methods  in  particular  are  very  flexi¬ 
ble  and  powerful  for  transport  calculations  in  this  semiclassical  regime.  Unfortunately,  these 
methods  fail  in  the  high-field,  high-scattering  rate  environments  which  are  of  increasing  interest 
in  microstructured  semiconductor  devices.54  In  situations  where  the  quantum  mechanical  nature 
of  the  system  becomes  important,  a  technique  which  has  the  power  to  calculate  electron  tran¬ 
sport  has  not.  yet  been  found.  Although  a  reasonably  complete  theoretic  framework  has  been 
establish  "  quantum  transport  has,  to  this  point,  proven  to  be  computationally  intractable. 

Collision  broadening  is  the  spread  of  the  electron  energy-momentum  relation  due  to  scatter¬ 
ing.  This  quasiparticle  nature  of  the  electron  is  best  described  by  the  finite  width  of  the  spectral 
function58  of  an  electron  interacting  with  its  environment  In  a  seminal  paper,  Chang  et  al.59 
explored  the  idea  of  including  collision  broadening  into  Monte  Carlo  calculations  of  electron 
transport  in  semiconductors.  Their  approach  was  to  include  the  broadening  of  the  electronic 
states  using  the  ideas  developed  by  Barker56  in  his  work  on  quantum  transport  and  the  mtracol- 
lisional  field  effect.  This  paper  by  Chang  et  al.59  has  received  much  comment  and  has  spurred 
further  research  into  the  problem.80,81  Recently,  Lugli  et  al.60  included  the  broadening  of  the 
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electronic  states  in  a  model  Monte  Carlo  calculation  using  a  more  precise  method,  which  makes 
use  of  the  many-body  spectral  function.  Their  work  used  a  single  parabolic  band  with  no  upper 
band  edge  and  the  Born  approximation  for  scattering. 

In  this  section,  a  method  which  includes  one  particular  quantum  effect,  collision  broaden¬ 
ing,  into  the  semiclassical  Monte  Carlo  method  by  applying  the  spectral  function  technique  is 
described  in  detail.  Comparisons  will  be  made  with  the  treatment  of  Chang  et  al.59  and  Lugli  et 
al.60  It  should  be  emphasized  that  this  treatment  does  not  constitute  a  complete  quantum  tran¬ 
sport  method.  It  incorporates  only  one  of  many  quantum  phenomena  into  a  semiclassical  tran¬ 
sport  calculation.  A  full  quantum  treatment  is  still  a  problem  for  the  future. 

The  method  discussed  here  is  based  on  the  semiconductor  model  that  was  employed  by 
Chang  et  al.59  The  full  band  structure  of  the  bulk  GaAs  conduction  band  is  calculated  using  an 
empirical  nonlocal  pseudopotential  method.62  The  band  structure  is  needed  for  a  quantitative 
and  even  qualitative  understanding  of  high-field  transport  in  semiconductors  because  of  the  large 
effect  of  the  density  of  states  on  the  scattering  rates.  Our  model  still  simplifies  the  electron 
dynamics  by  approximating  the  full  electron  scattering  at  high-electron  energies  from  all  the 
scattering  mechanisms  by  an  interaction  with  a  single  nonpolar  optical  phonon  mode  with  fre¬ 
quency  o;0.  The  coupling  between  the  electrons  and  the  phonons  is  chosen  so  that  the  cross  sec¬ 
tion  at  an  energy  of  0.6  eV  above  the  conduction  band  edge  is  equal  to  the  total  cross  section  of 
electrons  scattering  from  all  the  scatterers  in  GaAs.83  Thus,  the  well-known  analytic  formulas23 
are  used  up  to  0.6  eV  and  the  method  described  below  is  used  above  0.6  eV  for  the  calculation  of 
scattering  rates. 

The  electron  scattering  rate  due  to  the  phonons  is  calculated  in  the  Fock  approximation64 
by  self-consistently  solving  the  equation, 
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^  (2jt)3  E-^o-f(lT-P  f-E  (F-P  .E-lr&Jb)  (2'18) 

for  the  electronic  self-energy  S  which  gives  the  perturbation  of  electronic  states.  The  scattering 
rate  can  be  obtained  from  the  relation 

|--f  Im  S(E(%  (2.19) 

In  Eq.  (2.18),  g(F )  is  the  coupling  between  the  electrons  and  the  phonons,  and  e(T)  is  the  disper¬ 
sion  relation  for  the  electrons  as  given  by  the  band  structure  calculation.  The  use  of  Eq.  (2.18) 
for  the  self-energy  includes  the  approximation  that  phonon  absorption  processes  are  negligible. 
At  the  temperature  of  77  K  which  was  used  in  this  calculation,  the  phonon  occupation  coefficient 
for  absorption  is  of  the  order  10-2  while  the  emission  coefficient  is  of  order  1.  In  the  case  of  non¬ 
polar  optical  phonons,  g(T)  is  well  approximated  by  a  constant.  As  pointed  out  by  Chang  et 
al.,59  when  the  electron-phonon  coupling  is  independent  of  the  wavevector,  the  self-energy  is  a 
function  of  the  electron  energy  alone,  S(T,E)  =  E(E). 

The  electron  scattering  rates  calculated  by  the  conventional  analytic  method  and  the 
method  discussed  above  are  shown  in  Fig.  2.7  as  a  function  of  energy.  The  phonon  energy 
cu0  =  29  meV  and  the  electron-phonon  coupling  g  =  0.2  eV  have  been  used.  Due  to  the  decrease 
in  the  density  of  states,  the  scattering  rate  by  the  self-energy  method  starts  to  deviate  from  the 
rate  by  the  conventional  method  at  high  energies.  Figure  2.8  shows  the  self-energy  as  calculated 
self-consistently  from  Eq.  (2.18)  and  the  self-energy  calculated  by  the  Born  approximation  (where 
the  £  on  the  right-hand  side  of  Eq.  (2.18)  is  replaced  by  an  imaginary  infinitesimal  t<5)  In  the 
Born  approximation,  the  imaginary  part  of  the  self-energy  Im  £(E)  is  proportional  to  p0(E— fToj0) 
where  p0  is  the  density  of  states.  Thus  one  can  see  from  Fig.  2.8  that  the  conduction  band 
extends  to  approximately  4  eV.  It  is  also  evident  that  the  Fock  approximation  results  in  consid¬ 
erably  reduced  scattering  near  the  center  of  the  band  where  the  density  of  states  is  high.  This 


Figure  2.7 


Electron  scattering  rates  at  77  K  as  a  function  of  energy.  The  solid  line  is 
obtained  by  the  conventional  method  using  the  analytic  formulas,'3  and  the 
dashed  line  by  the  self-energy  method. 
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approximation  is  therefore  needed  for  a  more  accurate  calculation  of  the  scattering  rates  at  high 
energies. 

Collision  broadening  is  introduced  into  Monte  Carlo  calculations  by  the  inclusion  of  a  dis¬ 
tribution  for  the  electron  energy  after  each  collision.  In  the  usual  semiclassical  transport  calcula¬ 
tion,  the  rate  of  transition  of  the  electron  state  from  before  to  after  a  phonon  collision  is  given 
by  the  golden  rule: 


WfTTf.lTj)  ~  |V  |26(€(Trf)-c(iri)-h-a;0) 


(2.20) 


where  V  is  the  effective  electron-phonon  potential.  The  final  state  is  given  exactly  by  the  initial 
state  and  the  phonon  energy.  As  described  in  some  detail  by  Reggiani  et  al.,65  the  quasiparticle 
nature  of  the  interacting  electron  can  be  included  by  replacing  the  delta  function  in  the  golden 
rule  by  a  joint  spectral  density  K.  The  joint  spectral  density  K  is  a  convolution  over  the  initial 
and  final  state  spectral  functions, 


KflWT,)  =  /-^A(iTf,E-fra;0)A(iri,E). 


The  transition  rate  in  this  case  is  given  by 


(2.21) 


wfM)  ~  |v  |2K(iTf,ri). 


(2.22) 


In  Eq.  (2.21),  the  spectral  function  for  the  interacting  electron  is  calculated  from  the  self-energy 


a(Xe)  = _ -2  ImSpT.E) _ 

[E— €(1T)— ReE(ir,E)]2  +  (ImE^E))2 


(2.23) 


The  final  electron  state  after  scattering  by  a  phonon  is  given  by  the  joint  spectral  density  distri¬ 
bution  K(Trr,rs)  reflecting  the  spread  in  the  quasiparticle  energy.  It  is  obvious  that  this  method 
of  including  the  broadening  of  the  electronic  states  reduces  to  the  golden  rule  when  the  bare  elec¬ 
tron  spectral  function,  A(TT,E)  =  2?r<5(E — c(Tc)),  is  used  in  Eqs.  (2.21)  and  (2.22). 
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Since  the  self-energy  is  a  function  of  the  electron  energy  alone,  K  is  a  function  of  only  the 
initial  and  final  electron  energies,  K^k,)  =  K(f(lTf),£(lcj)).  The  quasiparticle  energies  E  are 
integrated  out  in  Eq.  (2.21).  In  Fig.  2.9,  the  joint  spectral  density  K(£f,£j)  is  plotted  as  a  func¬ 
tion  of  the  final  energy  £f  for  three  values  of  the  initial  energy  et.  As  can  be  seen  in  the  figure,  K 
is  asymmetric  with  an  enhanced  tail  which  extends  into  the  region  of  the  band  where  the  density 
of  states  is  largest.  This  tail  is  due  to  the  peak  in  the  density  of  states  near  the  band  center. 
For  £j  below  the  peak  in  the  density  of  states,  this  tail  in  K  extends  out  for  ef  >  £i;  while  for  £| 
above  the  peak,  this  tail  extends  out  for  £f  <  £;. 

It  is  also  important  to  note  that  the  integral  fd£ fK(€f,€j)  =  1.  This,  of  course,  is  needed  for 
conservation  of  the  electrons.  This  integral  over  K  gives  a  check  on  the  numerical  calculation. 
Finally  it  should  be  pointed  out  that  this  model  includes  the  real  part  of  E(E),  unlike  the  model 
of  Lugli  et  al.,80  which  gives  the  renormalization  of  the  electronic-band  energies  due  to  the  pho¬ 
nons.  This  makes  the  model  completely  consistent  with  the  quasiparticle  picture  of  the  elec¬ 
trons. 

To  study  the  importance  of  collision  broadening,  Monte  Carlo  calculations  of  the  electron 
energy  distribution  in  a  uniform  electric  field  have  been  performed  iD  GaAs  using  the  model 
described  above  at  77  K.  Polar  optical  phonons  and  intervalley  scattering  have  been  included 
for  energy  loss  of  low-energy  electrons  and  intervalley  electronic  motion  respectively.  The  possi¬ 
bility  of  impact  ionization  has  also  been  considered.  Figure  2.10  shows  the  steady-state  electron 
distribution  for  an  external  field  of  500  kV/cm.  Included  in  this  figure  are  results  of  calculations 
using  three  different  forms  for  the  transition  rate.  The  solid  curve  in  Fig.  2  10  gives  the  result 
for  a  simulation  using  the  golden  rule  transition  rate  as  given  by  Eq.  (2.20).  The  dotUci  curve 
gives  the  electron  distribution  with  the  full  collision  broadening  as  described  above  included. 
The  dashed  curve  is  a  simulation  using  a  Lorentzian  approximation  to  the  broadening  similar  to 
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that  used  by  Chang  et  al.59  In  this  approximation,  the  width  of  the  state  before  the  scattering 
ImE^i)  has  been  used  to  broaden  the  transition  rate.  For  this  initial  state  broadening,  the 
spectral  density  distribution  is 


K(€„€j)  = 


—2  ImE(€j) 


(2.24) 


This  approximation  simplifies  the  calculation  of  the  joint  spectral  density  K( ef,€j)  considerably 
while  including  an  approximate  broadening  of  the  electronic  states.  This  Lorentzian  state 
broadening  reproduces  many  of  the  features  of  the  full  broadening  calculation.  It  should  be  pos¬ 
sible  to  improve  on  this  approximation  by  including  the  spread  in  both  the  initial  and  final 
states  as  suggested  by  the  work  of  Barker.56 

It  is  easily  seen  from  Fig.  2.10  that  the  collision  broadening  has  a  profound  effect  on  the 
calculated  electronic  distribution.  Although  the  average  energies  of  all  three  models  are  equal 
within  the  accuracy  of  the  calculation,  the  simulations  which  included  broadening  have  a  much 
larger  spread  in  the  electronic  distribution  with  a  maximum  at  a  somewhat  lower  energy  and  a 
particularly  pronounced  enhancement  in  the  high-energy  tail.  This  is  because  with  the  collision 
broadening,  the  electrons  at  the  energy  range  where  the  density  of  states  is  high  (here,  between  1 
eV  and  2  eV)  tend  to  be  scattered  away  from  the  original  energies  due  to  the  large  spread  m  the 
spectral  density.  Thus,  the  broadened  distributions  have  a  very  nonthermal  (nonexponenual) 
high-energy  tail. 

Unlike  this  result,  Lugli  et  al.60  found  a  run-away  phenomenon  when  collision  broadening 
was  included.  The  average  electron  energy  of  their  broadened  distribution  was  significantly 
higher  than  that  calculated  using  the  golden  rule.  The  inclusion  of  the  precise  band  structure 
inhibits  this  run-away.  The  reason  for  this  can  be  seen  by  considering  the  function  K  plotted  in 
Fig.  2.9.  High  in  the  band  the  scattering  probability  is  skewed  towards  lower  energies  due  to 
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the  larger  density  of  states  there.  For  the  strictly  parabolic  dispersion  relation  used  by  Lugli  et 
al.,80  the  density  of  states  always  increases.  In  out  results,  the  electron  energy  distribution  is 
larger  in  both  the  high-  and  low-energy  regions,  than  the  golden  rule  results  with  little  or  no 
change  in  the  average  energy.  There  is  no  run-away  when  the  full  band  structure  of  the  semi¬ 
conductor  is  considered  in  the  electronic  model. 

The  high-energy  electron  population  enhancement  due  to  broadening  is  in  the  region  of  the 
band  where  processes  such  as  impact  ionization  begin  to  occur.  This  significant  change  in  the 
electronic  distribution  will  have  a  large  effect  on  the  calculated  rates  of  any  high-energy  effects. 
The  accurate  calculation  of  these  processes  would  seem  to  require,  therefore,  the  inclusion  of  the 
broadening  of  the  electronic  states  However,  caution  is  needed  for  the  numerical  application  of 
this  spectral  density  function  technique.  Since  the  spectral  density  is  not  a  classical  probability 
density,  one  often  has  to  make  approximations,  such  as  limiting  the  final  energy  within  certain 
energy  ranges  or  making  the  spectral  density  always  positive  definite,66  in  order  to  implement 
this  quantum  mechanical  quantity  into  a  semiclassical  calculation.  It  is  important  to  make  sure 
that  the  results  are  not  significantly  affected  by  these  approximations. 
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CHAPTER  3 

STUDY  OF  ELECTRON  IMPACT  IONIZATION  RATE  IN  GaAs 

3.1  Introduction 

Carrier  multiplication  by  impact  ionization  is  an  essential  mechanism  in  semiconductor 
devices  such  as  photodiodes  and  IMP  ATT  diodes,66  and  has  been  subject  to  many  theoretical 
and  experimental  investigations.  ’  ~  Most  of  the  analytical  theories  1  contain  several  adju¬ 
stable  parameters  whose  physical  meaning  are  not  easy  to  understand.  Recently,  a  Monte  Carlo 
algorithm  which  does  not  contain  adjustable  parameters  and  includes  a  realistic  band  structure 
was  developed  and  applied  by  Shichijo  and  Hess63  to  investigate  high  field  transport  in  GaAs. 
More  refined  Monte  Carlo  calculations  of  the  impact  ionization  rates  over  a  wide  range  of  electric 
fields  have  been  carried  out  by  Brennan  and  Hess.72  The  data  of  these  calculations63,72  are  in 
good  agreement  with  experimental  data,  for  example,  given  by  Bulman  et  al.75  In  these  calcula¬ 
tions,  a  steady-state  model  for  uniform  electric  fields  has  been  used  The  steady-state  model  is 
appropriate  for  analyzing  the  data  measured  in  the  p-i-n  diodes  or  punch-through  p-n  diodes 
with  long  depletion  regions.  However,  most  of  the  rate  measurements  have  been  carried  out 
either  in  abrupt  p+-n  junctions74,75  or  Schottky  barriers73  where  the  electric  fields  are  not  con¬ 
stant.  Also,  the  random  distribution  of  dopants  results  in  fluctuations  of  the  electric  fields  and 
affects  electron  transport.68  It  is,  therefore,  desirable  to  study  the  impact  ionization  in  nonuni¬ 
form  electric  fields,  which  can  includes  nonlocal  effects,  as  has  been  pointed  out  by  Capasso  et 
al76 

In  this  chapter,  the  electron  impact  ionization  in  GaAs  is  studied  in  various  situations  with 
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an  empinca1  oseudopotential  band  structure.  Included  scattering  mechanisms  are  polar  optical 
phonon  scattering  and  equivalent  intervalley  scattering  as  well  as  nonequivalent  intervallev 
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scattering.  Section  3.2  deals  with  the  dependence  of  the  impact  ionization  rate  on  the  details  of 
the  band  structure.  In  Section  3.3,  the  spatial  evolution  of  the  ionization  rate  and  the  average 
electron  energy  will  be  discussed  for  field  variations  typical  for  p+-n  junctions  or  Schottky  bar¬ 
riers.  The  effects  of  field  fluctuation  due  to  the  random  distribution  of  dopants  will  be  studied 
in  Section  3.4.  The  final  section  is  devoted  to  the  study  of  a  new  type  of  device  (avalanche  elec¬ 
tron  emitting  diodes77,78)  which  is  heavily  dependent  on  the  impact  ionization  process. 
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3.2  Band-Structure  Dependence 

In  the  previous  studies,83,72  the  empirical  local  pseudopotential  method  (as  developed  by 
Cohen  and  Bergstresser79)  was  used  to  obtain  the  band  structure  This  band-structure  calcula¬ 
tion  can  be  improved  by  using  the  empirical  nonlocal  pseudopotential  method  of  Chelikowsky 
and  Cohen62  in  which  the  effect  of  nonlocality  or  energy  dependence  of  the  pseudopotential  is 
included.  This  difference  in  method  brings  changes  in  the  band  structure,  and  affects  the  calcu¬ 
lation  of  impact  ionization  rates.  By  carefully  studying  these  changes,  one  can  get  important 
information  on  the  basic  characteristics  of  impact  ionization. 

Figure  3.1  shows  the  first  conduction  band  of  GaAs  calculated  by  the  local  and  nonlocal 
pseudopotential  method.  Compared  with  the  local  band  structure,  the  nonlocal  band  structure 
has  slightly  lower  energies  for  all  values  of  k.  This  difference  in  energy  between  the  two  bands  is 
more  pronounced  at  high  energies.  A  comparison  of  the  density  of  states  computed  with  both 
methods  is  shown  in  Fig.  3.2.  As  can  be  seen,  the  profile  of  the  nonlocal  method  exhibits  a 
higher  peak  in  the  density  of  states.  The  states  are  also  shifted  to  lower  energies  in  the  nonlocal 
calculation.  Since  the  scattering  rate  by  phonons  is  proportional  to  the  final  density  of  states, 
one  can  expect  that  the  energy  loss  of  electrons  will  be  higher  for  the  nonlocal  case  and  the  aver¬ 
age  energy  will  be  lower.  Accordingly  the  impact  ionization  rate  should  be  smaller. 
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Figure  3.2  Total  density  of  states  of  the  first  conduction  band  of  GaAs  as  a  function  of  electron 
energy  The  data  of  local  are  obtained  from  the  local  band  structure  in  Fig.  3.1,  while  the 
data  of  nonlocal  are  obtained  from  the  nonlocal  band  structure  in  Fig.  3.1. 


The  impact  ionization  rates  calculated  by  a  Monte  Carlo  method  with  these  two  band 
structures  are  shown  in  Fig.  3.3.  The  data  for  *h?  local  method  are  directly  from  Brennan  and 
Hess.72  As  expected  the  impact  ionization  rate  is  considerably  lower  for  the  nonlocal  method 
than  for  the  local  method  in  the  entire  range  of  electric  fields.  This  difference  is  relatively  small 
at  high  fields  and  increases  as  the  field  decreases.  The  reason  for  this  is  that  the  high-energy  tail 
of  the  distribution  is  sensitively  affected  by  the  details  of  the  band  structure  and  especially  by 
the  increase  in  the  density  of  states.  At  low  electric  fields,  it  is  this  high-energy  tail  (the  lucky 
and  nearly  lucky  electrons)  which  contributes  most  to  the  ionization  rate.  At  high  fields,  the 
distribution  function  is  heated  so  much  that  the  bulk  of  the  spherically  symmetric  part  of  the 
energy  distribution  function  contributes  to  the  ionization,  and  the  band  structure  has  less 
influence. 

Note  also  that  the  increased  scattering  for  the  nonlocal  band  structure  leads  to  the  disap¬ 
pearance  of  the  small  anisotropy  of  a  which  is  observed  for  the  local  band  structure.  In  other 
words,  increased  intervalley  scattering  shifts  the  anisotropy  to  much  lower  fields  for  the  nonlocal 
case. 

3.3  Ionization  Rates  in  Nonuniform  Fields 

There  are  two  widely  used  ways  to  define  the  impact  ionization  rate.  One  is  based  on  the 
rate  of  the  change  of  the  particle  density,  and  the  other  on  the  rate  of  the  change  of  the  current 
density.  Beni  and  Capasso80  have  shown  that  these  two  definitions  are  in  general  not  the  same 
and  the  difference  is  proportional  to  the  derivative  of  the  drift  velocity.  Of  course,  in  homogene¬ 
ous  steady-state  problems,  these  two  definitions  give  the  same  result.  However,  in  the  spatially 
transient  case,  the  difference  may  be  significant.  It  is  generally  accepted  that  the  definition 
based  on  the  current  density  must  be  used  in  the  transient  cases.80  Therefore,  the  connection 
between  these  two  definitions  and  the  microscopic  model  used  in  an  actual  Monte  Carlo 
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simulation  (the  inverse  mean  free  path  between  ionizing  collisions)  should  be  investigated  for  the 
appropriate  evaluation  of  the  impact  ionization  rate  in  nonlocal  problems. 

Thornber81  has  derived  the  ionization  rate  from  the  Boltzmann  equation,  and  has  shown  its 
equivalence  to  the  definition  of  the  ionization  rate  in  terms  of  the  current  density.  In  his  deriva¬ 
tion,  the  electron  impact  ionization  rate  is  give  by 


eehe 


=  S  /d3Pi/d3pj/d3pk/d3pi  RijkiCPi/Pj.PkiPi)  fi(%Pi,t)  /  Je&t) 
ijki 


(3.1) 


where  Je(7,t)  is  the  electron  particle  flux  density  and  Rjjki(l?j)lSj,]?k;]Ji)  is  the  rate  at  which  a  (l,p>j) 
carrier  scatters  to  (j^),  creating  an  electron-hole  pair  (i,"i?i)— (k,!^).  This  equation  can  be  written 


as 


(3.2) 


<*,(*,  t)  =  -7®r  E/d3Pi  fi(^,?i,t)  £/d3Pi/d3pj/d3pk  RijkiC?i,?j(^k;?i) 
Je(x,t)  1  ijk 

-  — ~r  S/d3 Pi  R.m)  lta,t) 

v„(3f,t)  i 


where  vn(?,t)  is  the  electron  drift  velocity.  The  function  f(7,{?i,t)  is  the  probability  distribution 
of  the  conduction  band  and  R|(i?i)  is  the  rate  at  which  a  (l,l>i)  electron  ionizes.  The  rate  per  unit 
time  of  creating  a  pair  by  impact  ionization  is  defined  as 


*■„(*, t) 


=  S/d3Pi  RiCft)  fOCpi.t)- 


(3.3) 


Then  Eq.  (3.2)  becomes 


«e(X.t)  = 


1 
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vnC^,t)rn(^.t)  Xn(3?,t) 


(3.4) 


which  is  precisely  the  inverse  mean  free  path  between  ionizing  collisions.  Hence,  the  microscopic 
model  used  in  Monte  Carlo  calculations,  Thornber’s  definition,  and  therefore  the  definition  of  the 
impact  ionization  rate  in  terms  of  the  current  density  are  all  equivalent.  This  means  that  we 


can  directly  compare  the  measured  data  with  the  calculated  results  by  the  Monte  Carlo  simula¬ 
tion.  Thornber’s  definition  (and  therefore  the  microscopic  model)  is  valid  as  long  as  the  underly¬ 
ing  Boltzmann  equation  picture  is  correct  and  applies  to  nonlocal  problems.  More  detailed  dis¬ 
cussions  on  the  applicability  of  this  model  can  be  found  in  Ref.  81. 

For  the  microscopic  model,  it  is  essential  to  know  Rijki(f?i,f?j,7?it;I?i)  and  f(T,f?j,t)  precisely. 
The  Keldysh  formula63  has  been  used  to  compute  R|(l?i).  However,  the  use  of  the  Keldysh  for¬ 
mula  makes  it  difficult  to  find  the  precise  final  states  of  electrons  and  holes.  Approximations 
used  to  restart  these  particles  at  zero  energy  and  neglect  created  (or  secondary)  electrons  are 
good  for  steady-state  situations.  This  is  because  the  actual  energies  of  final  states  are  close  to 
zero  in  general,  and  the  nonlocal  effects  of  secondary  electrons  are  not  in  the  steady-state  calcu¬ 
lations  by  definition  (device  length  is  assumed  to  be  infinite).  However,  for  the  transient  case, 
the  created  electrons  which  may  have  an  energy  distribution  quite  different  from  the  launched 
electrons  have  nonnegligible  effects.  It  is  then  necessary  to  keep  track  of  the  secondary  electrons 
to  obtain  the  exact  nonlocal  effects.  The  importance  of  secondary  electrons  becomes  significant 
when  the  impact  ionization  rate  is  high. 

The  electric  field  profiles  considered  in  this  section  are  shown  in  the  inset  of  Fig.  3.4. 
Under  the  depletion  approximation,  these  profiles  describe  the  electric  fields  in  abrupt  p+-n  junc¬ 
tions  or  Schottky  barriers  which  are  commonly  used  to  measure  the  impact  ionization  rate.  The 
zero  of  the  coordinate  system  x=0  is  located  at  the  junction  interface  at  which  the  Monte  Carlo 
electrons  are  launched.  Electrons  are  launched  with  zero  energy.  In  the  calculation,  the  empiri¬ 
cal  nonlocal  pseudopotential  band  structure  shown  in  the  previous  section  is  used.  The  material 
parameters,  such  as  the  ionization  threshold  and  coupling  constant  for  phonon  scatterings,  are 
the  same  as  in  Ref.  72.  Due  to  the  very  high  electric  field  and  low  doping  density  (~1016cm~3) 
in  the  n  region,  ionized  impurity  scattering  is  not  included  27 


Impact  Ionization  Rate  (1/cm) 


Distance  (A) 


Figure  3.4  Electron  impact  ionization  rates  in  <100>  direction  with  the  corresponding  field 
profiles  in  the  inset.  The  length  of  the  dead  space  is  approximately  700  A  for  (a) 
and  (b),  and  500  A  for  (c)  and  (d). 
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Figure  3.4  shows  the  electron  impact  ionization  rate  in  GaAs  with  the  corresponding  field 
profiles.  As  can  be  seen,  the  ionization  rate  is  very  low  for  a  considerable  distance  (the  dead 
space)  before  it  reaches  steady  state.  Similar  characteristics  are  observed  for  the  average  energy 
as  shown  in  Fig.  3.5.  Since  the  electrons  are  launched  at  zero  energy,  a  certain  distance  is 
needed  for  the  electron  gas  to  gain  enough  energy  to  impact  ionize.  The  dead  space  is  one  of  the 
typical  nonlocal  effects.  Grant82  expected  this  distance  to  be  the  distance  needed  for  electrons  to 
be  ballistically  heated  to  the  threshold  energy.  Since,  however,  electrons  lose  energy  through  ine¬ 
lastic  collisions  while  they  gain  energy  from  the  electric  field,  longer  distances  than  those 
obtained  from  the  simple  ballistic  transport  theory  are  computed.  A  precise  definition  of  the 
dead  space  is  not  straightforward  in  general.  Here,  we  define  it  as  the  distance  needed  to  reach 
the  maximum  of  the  impact  ionization  rate.  Since  the  transient  rates  in  the  steady  regime  after 
the  dead  space  are  in  good  agreement  with  the  steady-state  rates,  the  distance  as  defined  above 
is  equivalent  to  the  distance  needed  to  reach  the  steady-state  ionization  rate.  As  can  be  seen 
from  Fig.  3.4,  the  length  of  the  dead  space  decreases  as  the  maximum  field  increases  since  the 
electrons  are  heated  up  more  quickly.  This  length  corresponds  to  the  distance  which  is  needed 
to  ballistically  gain  about  3.4  eV  (which  then  corresponds  to  approximately  twice  the  impact 
ionization  threshold  energy  and  two  and  one  half  times  the  band  gap  energy  for  electrons  in 
GaAs).  The  rate  measured  in  p+-n  junctions  are  usually  plotted  vs.  1/E  with  a  dead  space 
correction83  estimated  by  Grant’s  method. 

Figure  3.6  shows  the  impact  ionization  rates  for  electric  field  profile  (a)  of  Fig.  3.4  pointing 
in  two  different  crystallographic  directions  (<100>  and  <  1 1 1  > ).  For  comparison,  the  steady- 
state  ionization  rates  are  included.  An  anisotropy  seems  to  develop  at  low  fields.  The  agree¬ 
ment  of  transient  and  steady-state  data  after  the  dead  space  means  that  the  simple  steady-state 
method  with  a  uniform  field  can  be  used  to  calculate  the  ionization  rate  in  this  range 
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Figure  3.5 


Average  electron  energy  in  <100>  direction.  The  electric  field  profile  of  each 
case  is  given  in  Fig.  3.4. 
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In  these  calculations,  the  electrons  have  been  launched  at  the  maximum  field.  In  p+-n  junc¬ 
tions,  the  field  does  not  increase  infinitely  fast  to  its  maximum  strength  because  of  the  small 
depletion  region  on  the  p+  side.  A  similar  slower  increase  will  occur  in  a  Schottky  barrier 
because  of  image  force  effects.  If  a  factor  of  100  is  assumed  for  the  ratio  between  doping  levels 
Na/Nd  (which  is  typical  for  p+-n  diodes),  an  additional  100  A  is  added  to  the  depletion  region 
for  profile  (a)  of  Fig.  3.4.  Figure  3.7  shows  the  result  including  this  correction.  Electrons  are 
launched  at  x=— 100  A  and  the  appropriate  field  profile  is  included  between  x=— 100  A  and  x=0 
along  with  profile  (a)  of  Fig.  3.4  for  positive  x  values.  As  expected,  the  ionization  rate  is  higher 
at  a  given  x  in  the  dead  space  and  reaches  the  peak  faster  for  this  case.  In  the  steady  regime, 
the  rates  of  two  cases  are  in  good  agreement. 

The  correction  by  including  secondary  electrons  has  also  been  estimated,  and  is  shown  in 
Fig.  3.7.  Since  the  ionizing  and  secondary  electrons  are  indistinguishable  after  the  ionization 
event,  one  need  only  follow  a  single  electron  after  the  event,  with  proper  weighting,  to  model  the 
entire  ensemble.  In  this  calculation,  the  excess  energy  of  an  ionizing  electron  is  randomly  distri¬ 
buted  between  possible  final  states  by  the  energy  conservation  law,  and  the  two  final  electrons 
are  assumed  to  have  identical  energy  and  momentum  for  simplicity.  Thus,  if  a  sample  electron 
experiences  an  ionizing  collision,  all  the  subsequent  contributions  by  that  electron  are  weighted 
by  the  factor  of  2.  This  correction  has  a  negligible  effect  on  the  dead  space  since  only  a  very 
small  number  of  secondary  electrons  are  created  in  this  region.  However,  as  the  number  of 
secondary  electrons  grows  by  reaching  the  regime  with  the  high  ionization  rate,  the  ionization 
rate  becomes  smaller  than  that  without  the  correction  because  the  secondary  electron  distribu¬ 
tion  is  cooler  than  the  launched  electrons  with  a  dead  space  for  each  secondary  electron.  To 
accurately  account  for  this  effect,  the  Keldysh  formula  must  be  replaced  by  a  more  precise 
theory  which  can  give  Rijki("pi,"Pj>lA;f>i)  rather  than  R|(pj),  as  suggested  by  Kane.84 
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Figure  3  7 


Effects  of  the  corrections.  The  do  correction  case  corresponds  to  profile  (a)  of  Fig 
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3.4  Effects  of  Random  Distribution  of  Dopants 

The  effects  of  field  fluctuations  due  to  the  random  distribution  of  dopants  on  impact  ioni¬ 
zation  were  first  addressed  by  Shockley  in  his  celebrated  paper  on  p-n  junctions.68  Shockley 
assumed  that  dopants  are  incorporated  in  the  junction  independently,  and  therefore  follow  the 
Poisson  distribution.  This  assumption  is  justified  because  dopants  are  usually  incorporated  dur¬ 
ing  growth  or  processing  under  the  condition  of  high  temperature.  At  these  temperatures,  the 
thermally  generated  intrinsic  electron  and  hole  concentrations  are  high  enough  to  screen  the  elec¬ 
trostatic  interactions  between  dopants.  Portraying  dopants  as  a  continuum  with  fluctuating 
density,  Shockley  calculated  the  root-mean-square  voltage  deviation  between  two  points 
separated  by  a  distance  d,68 

<5V  =  [-a^No+NJdl1^  (3.5) 


Here  q  is  the  magnitude  of  the  electronic  charge,  e  is  the  permittivity,  and  ND  and  NA  are  aver¬ 
age  donor  and  acceptor  concentrations.  Equation  3.5  is  derived  by  considering  an  infinite 
depleted  semiconductor  with  uniform  average  doping.  Since  only  the  high  field  region  of  the  p-n 
junction  is  of  interest  in  connection  with  impact  ionization  phenomena,  it  is  valid  for  Shockley 
to  neglect  the  effects  of  the  edge  of  the  depletion  region.  From  Eq.  (3.5),  Shockley  proceeds  to 
derive  an  effective  electric  field  for  impact  ionization, 


(3  6) 


Here  F  is  the  electric  field  and  Ff  is  a  quantity  representing  the  fluctuations  in  the  field  which  is 
given  by 


- 

Eth/q 


(3.7) 
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where  Eth  is  the  threshold  energy  for  impact  ionization.  In  deriving  Eqs.  (3.6)  and  (3.7),  Shock- 
ley  assumes  that  the  average  electric  field  is  constant,  which  is  a  rather  weak  approximation. 
Since  the  effective  field  is  larger  than  the  average  field,  Shockley  predicted  strong  enhancement 
in  the  impact  ionization  rate  due  to  the  field  fluctuations. 

More  accurate  results  can  be  obtained  with  the  use  of  a  computer.  With  the  location  of 
each  discrete  dopant  randomly  chosen  in  a  p-n  junction,  the  realistic  field  fluctuations  can  be 
calculated  by  solving  the  Poisson  equation  in  three  dimensions  using  a  method  such  as  the  finite 
difference  SOR  algorithm.85  Treating  the  impurities  as  infinitesimal  point  charges  would  result 
in  singularities  in  the  solution  of  the  Poisson  equation.  These  singularities  would  contribute  to 
the  calculated  fluctuations  (which  are  still  finite)  in  an  unphysicaJ  way  if  included  into  the  field 
fluctuations  without  modification,  since  the  effect  of  the  core  electrons  is  totally  neglected. 
Thus,  it  is  more  appropriate  to  use  a  pseudopotential  and  therefore  cut  off  the  potential  varia¬ 
tions  within  the  first  Bohr-radius.  For  simplicity  of  calculation  in  this  model,  it  is  assumed  that 
the  carriers  which  contribute  to  the  current  do  not  significantly  affect  the  potential  energy 
profile.  It  is  also  assumed  that  the  applied  bias  affects  only  the  average  potential  energy  profile, 
but  not  the  fluctuations. 

The  particular  case  considered  is  a  GaAs  junction  doped  with  3.6  X  1018  cm-3  acceptors  on 
one  side  and  3.6  X  1018  cm-3  donors  on  the  other  side.  The  root-mean-square  voltage  fluctua¬ 
tions  between  two  points  in  space  is  plotted  in  Fig  3.8  as  a  function  of  the  distance  between  two 
points.  Shockley’s  result,  given  by  Eq.  (3.5),  is  shown  also  for  comparison.  The  curves  show 
reasonably  close  agreement.  In  the  calculation,  the  spatial  coordinate  normal  to  the  p-n  junction 
has  been  denoted  by  x,  and  the  transverse  coordinates  by  y  and  z  The  normal  component,  of 
the  electric  field  Fx(x,y,z)  is  obtained  by  numerically  differentiating  the  potential  Throughout 
this  section,  the  normal  component  of  the  electric  field  will  be  referred  to  as  the  field  in  the 
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remainder  of  the  text.  The  fluctuations  in  the  field  are  given  by 


<5Fx(x,y,z)  =  Fx(x,y,z)  -  Fx(x). 


(3.8) 


Here  Fx(x)  is  the  field  averaged  over  y  and  z.  The  fluctuations  are,  then,  superimposed  on  the 
field,  averaged  over  y  and  z,  for  a  p+-n  junction  under  high  reverse  bias  in  order  to  calculate  the 
electron  impact  ionization  rates.  This  average  field  is  calculated  from  the  depletion  approxima¬ 
tion,  which  gives  excellent  agreement  with  more  exact  calculations  of  the  average  electric  field 
except  near  the  edge  of  the  depletion  region.  Since  the  latter  is  a  low-field  region,  it  is  not 
important  to  impact  ionization.  Statistical  analysis  is  performed  on  one-dimensional  cuts  of  the 
field  as  a  function  of  x,  with  y  and  z  held  constant.  The  probability  density  of  the  root-mean- 
square  (rms)  field  fluctuations  for  the  one-dimensional  segments  can  be  computed  from 


<5F(y,z) 


[  Fx(xj,y,z)  -  F X(xj)  ]2 
N— 1 


(3.9) 


From  Eq.  (3.9),  one  can  see  how  much  the  field  fluctuates  in  the  normal  direction  to  the  junc¬ 
tion. 


To  save  computer  time  the  electron  impact  ionization  rates  are  not  calculated  for  the  entire 
three-dimensional  volume.  Two  one-dimensional  electric  field  profiles  are  chosen  from  the  proba¬ 
bility  density  of  the  rms  fluctuations  of  Eq.  (3.9):  one  for  which  the  standard  deviation  is  greater 
than  50%  of  those  for  one-dimensional  cuts  in  the  sample  mesh,  and  one  for  which  the  standard 
deviation  is  greater  than  90%  The  former  represents  a  typical  example  and  the  latter  an 
extreme  case.  These  two  electric  field  profiles  are  shown  in  Fig  3.9  along  with  the  average  elec¬ 
tric  field  given  by  the  depletion  approximation  for  a  reverse  bias  of  24  volts 

The  impact  ionization  rates  in  GaAs  for  the  electric  fields  given  in  Fig  3.9  are  shown  in 
Fig.  3.10.  The  same  transport  model  and  the  same  empirical  band  structure  are  used  as  given  in 
Section  3.3.  The  effect  of  the  secondary  electrons  is  not  included  since  it  is  negligible  as 


Figure  3.9  Electric  field  profiles  in  a  p+-n  junction.  90%  represents  the  case  for  which  the 
standard  deviation  is  greater  than  90%  of  those  for  one-dimensional  cuts  in  the 
sample  mesh,  and  50%  represents  the  case  for  which  the  standard  deviation  is 
greater  than  50%. 
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Figure  3.10  Electron  impact  ionization  rates  in  GaAs  in  <100>  direction  The  electric  field 
profile  for  each  case  is  given  in  Fig.  3.9. 
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discussed  in  the  previous  section.  In  the  dead  space  where  the  ionization  rates  are  very  low  in 
spite  of  the  high  field,  the  ionization  rates  are  almost  the  same  in  the  three  cases  even  though 
the  local  fields  have  quite  different  values.  This  is  because  the  electron  distribution  is  too  cold  to 
respond  to  the  variation  in  the  local  field  sensitively  and,  therefore,  the  nonlocal  effects  are  dom¬ 
inant.  When  the  electrons  gain  sufficient  energy  by  reaching  the  steady  region  after  the  dead 
space,  the  ionization  rates  become  affected  by  the  local  fluctuations  resulting  in  a  spatially 
fluctuating  ionization  rate.  When  the  electric  field  varies  slowly  in  space  with  respect  to  the 
mean  free  path  of  an  electron,  the  transient  ionization  rates  depend  only  on  the  local  field  and 
not  explicitly  on  the  position.88  As  discussed  in  the  previous  section,  such  is  the  case  in  the 
depletion  approximation.  However,  the  field  fluctuations  shown  in  Fig.  3.9  are  not  slowly  vary¬ 
ing.  The  response  of  electrons  to  these  large  field  gradients  is  delayed.  Therefore,  adjacent  posi¬ 
tive  and  negative  field  fluctuations  tend  to  cancel.  Hence,  the  big  fluctuations  such  as  observed 
near  x  =  1500  A  in  the  90%  case  in  Fig.  3.9  result  in  a  delayed  and  relatively  small  change  in 
the  ionization  rate  near  x  =  2000  A  (Fig.  3.10). 

This  nonlocal  effect  is  clearly  shown  in  Fig.  3.11  for  the  average  electron  energy.  To  gain 
energy  an  electron  must  travel  a  distance  through  a  field  while  suffering  from  inelastic  collisions. 
Since  the  fluctuations  are  very  rapid  when  compared  to  the  ionization  mean  free  path,  the 
energy  gain  of  the  electrons  from  the  ups  and  downs  of  the  field  are  virtually  averaged  out  and 
the  average  electron  energy  in  Fig.  3.11  shows  rather  small  changes.  These  small  changes  in 
average  energy  also  result  in  small  fluctuations  of  the  ionization  rate  as  discussed  above.  There¬ 
fore,  a  significant  increase  in  the  ionization  rate  due  to  the  field  fluctuation  is  not  observed, 
unlike  in  Shockley’s  steady-state  analysis.68  This  is  different  from  the  case  of  impact  ionization 
in  a  superlattice  where  the  changes  in  potential  energy  are  large  and  vary  less  rapidly  compared 
to  the  mean  free  path.12  No  averaging-out  will  occur  if  electrons  gain  a  considerable  amount  of 
energy  within  a  very  short  time  period  (<  rph  where  rph  is  the  phonon  scattering  time  constant). 
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3.5  Avalanche  Electron  Emitting  Diodes 

The  recent  report  of  efficient  silicon  cold  cathodes  by  van  Gorkom  et  ai. 77,78  has  triggered 
renewed  interest  in  such  solid-state  electron  sources.  Current  densities  as  high  as  1000  A/cm2 
and  reduced  brightnesses  of  900  A/(cm2  sr  V)  for  a  1  /im  diameter  emitting  area  have  been 
reported,  making  these  devices  attractive  for  use  in  numerous  applications.  Although  p-n 
emitters  were  proposed  long  ago,  only  recently  have  they  been  proven  practical  with  high- 
current  density  and  efficiency.  The  device  consists  of  a  shallow  p+-n++  junction  fabricated  at 
the  semiconductor  surface  with  the  n++  channel  at  the  surface  as  shown  in  the  inset  of  Fig.  3.12. 
Cesium  is  typically  absorbed  onto  the  surface  in  order  to  lower  the  work  function  to  approxi¬ 
mately  1.7  eV.  Reverse  biased  in  avalanche  breakdown,  electrons  entering  the  depletion  region 
from  the  p-type  side  of  the  junction  gain  energy  from  the  field  and  impact  ionize,  creating  free 
electron-hole  pairs.  Due  to  the  very  high  electric  field,  the  avalanche  multiplication  (i.e.,  impact 
ionization)  is  an  essential  part  of  the  device  operation  in  this  structure.  Electrons  which  have 
sufficient  energy  upon  reaching  the  surface  may  then  overcome  the  work  function  and  be  emitted 
into  the  vacuum.  The  efficiency  is  defined  as  tj  =  jyac/jjunc  where  jvac  is  the  current  density  in 
the  vacuum  and  jjunc  is  the  current  density  through  the  junction,  normal  to  the  surface.  The 
maximum  field  in  the  junction  is  typically  1.5X106  V/cm.77 

The  study  of  these  devices  is  interesting  from  a  practical  point  of  view  in  an  effort  to 
design  more  efficient  devices,  and  also  from  the  point  of  view  of  high-field  transport  physics.  In 
this  section,  the  potential  performance  of  GaAs  devices  is  studied.  The  results  are  compared  to 
the  calculated87  and  the  experimental77  data  for  Si.  The  transport  model  used  in  these  calcula¬ 
tions  is  essentially  the  same  as  described  in  the  previous  sections.  Unlike  GaAs,  the  first  two 
conduction  bands  are  included  in  the  Si  model,87  since  the  second  conduction  band  is  degenerate 
with  the  first  at  the  X  point  and  can  therefore  have  a  significant  population  at  high  electric 
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fields.88  In  GaAs,  the  second  conduction  band  lies  at  least  1  eV  above  the  first  and  would  have 
a  correspondingly  small  population.  Even  though  the  devices  are  heavily  doped  on  both  the  n- 
and  p-type  sides  of  the  junction,  ionized  impurity  scattering  is  not  included.  Since  only  the 
high-energy  electrons  are  of  interest,  and  for  these  electrons  Coulombic  scattering  prefers  small 
angles,  the  momentum  distribution  will  not  be  disturbed.  It  is  assumed  that  the  effect  of  the 
holes  on  the  avalanche  process  is  to  simply  create  additional  electrons  near  the  p-type  side  of  the 
junction,  which  will  behave  identically  to  the  thermally  generated  electrons  leaving  the  p-type 
region.  Under  these  assumptions,  the  electrons  created  by  hole-initiated  ionization  are  accounted 
for  in  the  model.  Band-to-band  tunneling  has  not  been  included  in  the  model,  though  a  10% 
tunneling  contribution  has  been  estimated  using  an  analytic  model.77  In  this  section,  the  secon¬ 
dary  electrons  created  by  impact  ionization  events  are  accounted  for,  as  described  in  Section  3.3, 
to  include  the  effect  of  avalanche  multiplication  on  both  the  number  and  distribution  of  elec¬ 
trons  reaching  the  barrier  at  the  surface.  In  addition,  a  numerical  enhancement  of  the  high- 
energy  tail  of  the  distribution89  has  been  used  to  speed  convergence  of  the  calculation  of  r). 

Finally,  the  transmission  of  the  electrons  over  the  barrier  into  the  vacuum  must  be  treated. 
For  the  classical  case  (i.e. ,  parabolic  energy  bands),  one  can  define  a  directional  component  of  the 
energy  such  as  Ex  oc  kx2.  The  quantum  treatment  of  the  transmission  would  require  a  Bloch 
function  description  of  the  electrons  in  the  semiconductor41  with  the  added  complications  of  sur¬ 
face  effects.  Instead,  a  simple,  analytic  model  for  the  transmission  coefficient  has  been  used  in 
this  treatment.  In  order  to  account  for  the  fact  that  not  all  electrons  with  total  energy  E^  can 
surmount  the  barrier,  we  have  used  an  effective  transmission  probability 
A  f  <9Er 

T.fKEt)  =  1  “  -grp  %  >  <t>  ,  >  0  (3.10) 

where  <f>  is  the  work  function  and  k^  is  the  component  of  the  crystal  momentum  normal  to  the 


surface.  This  treatment  is  exact  in  the  limit  of  a  spherically  symmetric  distribution  function  and 
parabolic  energy  bands.  Qualitative  results  such  as  the  comparison  of  different  materials  and  cal¬ 
culated  trends  should  not  be  seriously  affected  by  assuming  Eq.  (3.10).  Comparison  with  experi¬ 
ment  shows  that  a  constant-factor  scaling  of  the  efficiency  is  sufficient  to  obtain  good  agreement 
with  experiment. 

In  Fig.  3.12,  results  of  the  Monte  Carlo  calculations  are  compared  with  published  experi¬ 
mental  data77  for  t)  vs.  <j>.  For  the  experimental  data,  the  best-fit  function  r?  =  0.6exp(— d>/0.41) 
(where  <t>  is  in  eV’s)  has  been  used.  For  GaAs  the  device  parameters  used  in  the  calculation  are 
Na  =  1X1018  cm-3  and  Nq  =  5X1018  cm-3  The  Si  devices  are  assumed  to  have 
Na  =  1X1018  cm-3  and  Nq  =  5X1019  cm-3.  In  each  case,  the  device  is  assumed  to  be  biased 
such  that  the  maximum  electric  field  is  1.5X106  V/cm.  Corresponding  depletion  widths  are 
Wp=1070  A  and  Wn=210  A  for  GaAs,  and  Wp=980  A  and  Wn=20  A  for  Si.  The  undepleted 
n++  channel  width  d  is  70  A  for  both  the  GaAs  and  Si  devices.  These  values  of  d  are  chosen 
based  on  the  silicon  emitters  which  have  been  demonstrated.  The  results  show  that  GaAs  has  a 
much  smaller  calculated  efficiency  than  Si.  The  calculated  data  for  Si  follow  a  curve  which  is 
nearly  parallel  to  the  experimental  data  The  difference  in  magnitude  may  be  partially  attri¬ 
buted  to  the  assumption  of  Eq.  (3.10)  for  the  transmission  coefficient.  However,  the  major  rea¬ 
son  for  the  difference  appears  to  be  of  experimental  origin:  it  is  known  that  the  measured 
efficiency  is  highly  sensitive  to  the  surface  contamination  In  a  more  recent  experiment,  emitters 
with  clean  surfaces  show  an  efficiency  of  4%~5%  and  up  to  8 %  at  <f>  =  1  7  eV.90  The  efficiency 
is  also  strongly  dependent  on  the  actual  dopant  distribution.  Diffusion  of  the  dopant  distribu¬ 
tion  results  in  lower  maximum  electric  field  and,  thus,  smaller  experimental  efficiencies.90  Con¬ 
sidering  these  two  points,  we  conclude  that  the  calculated  data  represent  an  upper  boundary 
value.  For  GaAs,  the  calculated  efficiency  is  also  nearly  parallel  to  the  experimental  Si  data.  We 
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would  expect  that  the  difference  between  the  measured  efficiencies  would  be  consistent  with  the 
difference  in  the  calculated  values,  so  the  Monte  Carlo  simulation  indicates  roughly  an  order  of 
magnitude  greater  efficiency  for  comparable  Si  sources. 

The  calculated  efficiency  is  plotted  as  a  function  of  d  for  both  GaAs  and  Si  in  Fig.  3.13. 
Note  that  for  both  materials,  the  efficiency  is,  to  a  good  approximation,  exponential  in  d.  The 
steeper  decrease  of  the  efficiency  for  GaAs  is  indicative  of  the  higher  scattering  rates  (larger  den¬ 
sity  of  states).  The  control  of  d  in  fabricating  the  device  is  therefore  even  more  important  for 
GaAs  devices.  Note  that  the  higher  n-type  doping  levels  possible  in  Si  allow  one  to  make  the 
channel  narrower  than  for  GaAs  while  maintaining  the  same  current-carrying  capability. 

For  both  Si  and  GaAs,  the  work  function  is  larger  than  the  threshold  energy  for  impact 
ionization,  <j>  >  Etj,.  Thus,  it  is  clear  that  impact  ionization  is  at  once  essential  for,  and  detri¬ 
mental  to,  the  operation  of  the  device.  The  junctions  must  be  biased  in  avalanche  breakdown 
for  a  sufficient  number  of  electrons  to  reach  the  surface  to  obtain  practical  values  of  jvac,  and  yet 
by  impact  ionizing,  many  electrons  lose  nearly  all  of  their  energy  and  are  then  unable  to  sur¬ 
mount  the  barrier.  A  soft  threshold  for  impact  ionization  is  necessary  to  enable  a  significant 
number  of  electrons  to  remain  in  the  high-energy  (Ep  >  <A)  tail  of  the  distribution.  From  various 
investigations,63,88,91  it  can  be  seen  that  the  threshold  in  Si  is  softer  than  in  GaAs.  Also,  the 
higher  density  of  states,  and  therefore  higher  scattering  rates,  at  energies  near  <±>  in  GaAs 
decrease  the  fraction  of  electrons  in  the  high-energy  tail.  These  effects  both  contribute  to  the 
lower  efficiency  of  GaAs  as  indicated  in  Figs.  3.12  and  3.13.  The  device  performance  could  be 
also  enhanced  by  relocating  the  generation  of  secondary  electrons  nearer  to  the  p— type  side  of 
the  junction.  Present  devices  have  not  been  optimized  in  this  respect. 
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CHAPTER  4 

ELECTRON  TRANSPORT  ACROSS  HETERO  JUNCTION  BAND  DISCONTINUITIES 

IN  STRONG  ELECTRIC  FIELDS 


4.1  Introduction 

I 

There  have  been  significant  theoretical  efforts  to  understand  the  effects  of  heterointerfaces 
on  the  carrier  transport.92-98  Most  of  the  efforts  have  been  devoted  to  the  transport  properties 
in  a  direction  parallel  to  the  interfaces  such  as  quasi-two-dimensional  transport  in  field  effect 
devices.92,93  However,  the  study  of  carrier  transport  across  heterointerfaces  is  also  essential  for 
the  understanding  of  many  interesting  nonlinear  effects  such  as  ionization  enhancement  in 
heterostructure  avalanche  photodiodes  (APD’s). 12,13 

In  this  chapter  a  detailed  study  of  nonlinear  electron  transport  across  heterointerfaces  in 
strong  electric  fields  is  presented,  based  on  an  ensemble  Monte  Carlo  simulation  method.  Unlike 
most  of  the  previous  studies  on  this  subject,94-96  particular  attention  has  been  devoted  to 
microscopic  features  such  as  the  evolution  of  the  energy  and  momentum  distribution  function  in 
space  and  time,  from  which  one  can  deduce  macroscopic  results.  It  is  shown  that  the  combined 
effects  of  real-space  and  momentum-space  structure  bring  pronounced  peaks  and  anisotropy  in 
energy  and  momentum  distribution  respectively.  The  dependence  of  the  ionization  rate  on  the 
band  structure  of  the  neighboring  layer  is  also  discussed. 


•1 


4.2  Description  of  Model 

Monte  Carlo  simulations  of  electron  transport  have  been  performed  for  various 
configurations  and  material  compositions  of  the  AlGaAs/GaAs  ?vstem.  F’igure  4  1  shows  two 
single-quantum-well  heterostructures  in  which  detailed  studies  have  been  carried  out.  The  struc¬ 
tures  shown  in  Fig.  4.1(a)  and  (b)  may  represent  part  of  the  i  region  of  a  reverse-biased 
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quantum-well  APD12  and  the  p+-i-n+-i  region  of  a  reverse-biased  doped  quantum-well  APD13 
respectively.  Generally  superlattice  structures  are  studied  in  connection  with  novel  versions  of 
these  devices.  However,  only  one  period  has  been  simulated  in  this  chapter  in  order  to  clarify 
the  fundamental  physics  of  nonlinear  transport. 

Each  layer  in  Fig.  4.1  is  assumed  to  be  ideal  and  free  of  defects  and  oriented  along  the 
<100>  direction.  Electrons  are  generated  in  a  region  prior  to  layer  1  by  an  external  source 
such  as  light  and  injected  into  the  layers  that  are  depleted  by  a  reverse  bias.  To  rule  out  the 
possible  nonlocal  effects  in  the  energy  distribution  when  electrons  enter  the  well  layer  from  layer 

o 

1,  a  length  of  1500  A  has  been  chosen  for  the  width  of  layer  1.  The  energy  and  momenta  of  the 
generated  electrons  are  assumed  to  obey  a  Maxwell-Boltzmann  distribution  before  they  are 
accelerated  by  the  high  fields.  Since  the  semiclassical  transport  is  intended  to  be  studied,  the 
well  layer  is  rather  long  (a  width  of  400  A  is  assumed  for  the  GaAs  layer  unless  otherwise 
specified).  Therefore,  size  quantization  and  resonance  effects  are  of  minor  importance.  The 
actual  size  of  layer  3  is  unimportant  for  the  results  which  we  report.  Unlike  in  Fig.  4.1(a),  the 
electric  field  in  Fig.  4.1(b)  is  not  uniform  The  electric  field  in  layer  1  varies  from  500  kV/cm  to 
100  kV/cm  due  to  the  change  in  the  build-in-field,  while  the  field  in  layers  2  and  3  is  fixed  to 
100  kV/cm.  Such  field  distributions  could  be  obtained  in  multiquantum-well  structures  by 
means  of  aD  appropriate  periodic  doping  profile.13  The  A1  mole  fraction  for  layers  1  and  3  is 
chosen  to  be  0.36  unless  otherwise  specified.  To  get  the  conduction  band  offset,  the  70/30  rule9' 
has  been  used.  The  empirical  band  structure  in  GaAs  is  calculated  as  in  Chapter  3.  while  for 
AlGaAs  the  "analytical"  band  structure  as  described  in  Ref.  98  is  used.  Since  the  electrons  in 
the  AlGaAs  layers  are  hot  due  to  the  high  electric  field,  inclusion  of  a  more  realistic  (but 
currently  not  available)  band  structure  for  AlGaAs  may  be  of  importance  The  material  param¬ 
eters  used  in  this  study  are  summarized  in  Table  4.1. 
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Table  4.1 


Material  parameters  for  the  simulation  program 


Effective  mass  ratio 
(m/m0) 

Nonparabolicity 

(eV-1) 

Valley  separation  relative 
to  the  T-valley  minimum 
(eV) 

• 

GaAs 

T  valley 

L  valley 

0.069 

1.47  (longitudinal) 

0.12  (transverse) 

0.67 

0.44 

0.29 

X  valley 

1  58  (longitudinal) 

0.24  (transverse) 

0.36 

0  48 

• 

F  valley 

0.069+0. 103*x 

0.67-0.94*x 

AlxGaj_xAs 

L  valley 

1,47  (longitudinal) 
0.12+0.036*x  (transverse) 

0.44-0.038*x 

0.29-0. 647  *x 

• 

X  valley 

l.o8+1.21*x  (longitudinal) 
0.24  (transverse) 

0.36 

0.48-1. 156*x 

Included  scattering  mechanisms  are  polar  optical  and  deformation  potential  scattering  as 
well  as  impact  ionization.  Electron-electron  scattering  is  also  included  in  the  well  layer  to 
account  for  the  charge  accumulation  which  can  occur.  Impact  ionization  rates  are  calculated  by 
using  the  estimator  proposed  by  Thornber81  as  discussed  in  Chapter  3.  The  ionization  threshold 
energy  and  the  p  factor  (of  the  Keldysh  formula)  for  GaAs  are  chosen  to  be  1.7  eV  and  0.5 
respectively.  For  AlGaAs,  the  threshold  energy  has  been  chosen  proportional  to  the  increase  in 
the  energy  gap  over  GaAs.  In  the  calculation,  the  effect  of  secondary  electrons  has  also  been 
included  for  some  extent.  For  simplicity,  the  ionizing  and  ionized  (secondary)  electrons  are  given 
zero  energy  after  ionization  occurs. 

At  the  interfaces,  electron  transfer  for  the  L  and  X  valley  electrons  is  treated  classically  by 
conserving  electron  energy  and  parallel  momentum.  With  classical,  it  is  meant  that  an  electron 
can  cross  the  interface  only  when  it  has  high  enough  energy  to  overcome  the  barrier.  This 
simplification  is  motivated  by  the  fact  that  the  time  constants  for  tunneling  are  relatively  long 
compared  to  the  short  transients  due  to  the  high  field.  For  the  T  valley  electrons,  the  tunneling 
process  is  included  along  with  the  classical  transfer.  By  conserving  energy  and  particle  flux,  the 
WKB  approximation  is  used  to  calculate  the  tunneling  probability.  The  reflections  at  the  inter¬ 
faces  are  assumed  to  be  specular.  The  quantum  interference  effects  by  the  reflections39  are  unim¬ 
portant  since  the  inelastic  mean  free  path  is  much  shorter  than  the  well  width. 

In  the  calculation,  a  constant  time  step  (2  X  10-15  sec)  discretization  scheme  is  used.  This 
allows  the  tracking  of  the  time  evolution  of  electron  transport  parameters.  For  the  electric  field 
profile,  the  self-consistent  correction  is  not  included.  At  high  fields  such  as  used  in  these  calcula¬ 
tions,  the  band  distortion  by  the  trapping  effects  in  the  well  is  not  important.12 
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4.3  Results  and  Discussion 

While  crossing  the  interface  and  propagating  into  the  GaAs  layer,  electrons  gain  excess 
energy  from  the  conduction  band  offset  and,  therefore,  exhibit  an  overheated  energy  distribution. 
F'gure  4.2  shows  the  electron  energy  distributions  with  an  electric  field  of  500  kV/cm  in  struc¬ 
ture  (a)  of  Fig.  4.1.  The  energy  distributions  in  layers  1  and  3  are  calculated  with  40  A  from 
the  interface  while  the  displayed  distribution  for  layer  2  represents  an  average  over  the  whole 
layer.  In  this  figure,  the  energy  distribution  in  layer  2  is  shifted  to  higher  energy  (by  the 
amount  of  the  band  edge  discontinuity)  as  expected,  and  exhibits  two  peaks,  one  at  very  low 
energy  near  the  T  valley  minimum  and  the  other  at  high  energy.  The  reason  for  the  low-energy 
peak  is  two-fold.  The  shift  in  energy  distribution  enhances  the  impact  ionization  and  creates 
cooled  electrons  in  the  T  valley.  Also,  the  heterointerface  barriers  tend  to  confine  (by  reflection) 
cold  electrons  in  the  well  region.  However,  most  of  the  cold  electrons  are  reheated  quickly  due  to 
the  high  held  applied  and,  thus,  the  low-energy  peak  is  reduced  in  spite  of  the  high  ionization 
rate  As  the  field  decreases,  the  low-energy  peak  also  decreases  along  with  the  ionization  rate  as 
shown  in  Fig.  4.3.  The  confinement  effect  is  negligible  in  this  case. 

The  overheating  causes  strongly  enhanced  electron-phonon  interactions.  Therefore,  elec¬ 
trons  in  layer  2  relax  to  the  steady-state  distribution  rather  rapidly.  As  shown  in  Fig.  4.4,  this 
quick  relaxation  results  in  nearly  the  same  distribution  as  in  the  situation  illustrated  in  Fig.  4.2. 
although  the  initial  distribution  in  layer  1  is  different  due  to  the  different  A1  mole  fraction.  The 
small  peaks  for  layers  1  and  3  are  due  to  the  low  energetic  electrons  which  overcome  or  tunnel 
through  the  potential  barrier 

For  the  simulation  of  structure  (b)  in  Fig.  4.1,  however,  electrons  experience  strong  relaxa¬ 
tion  and  confinement  due  to  the  relatively  small  electric  field  in  the  GaAs  layer  (100  kV/cm)  In 
this  case,  a  smooth  distribution  enters  the  GaAs  layer  from  layer  1  and  relaxes  to  a  distribution 


Distribution  (arbitrary  units) 


Distribution  (arbitrary  units) 


Figure 


4.4  Electron  energy  distributions  in  structure  (a)  of  Fig.  4.1.  The  electric  field  is  500 
kV,  cm.  The  A1  mole  fraction  is  0.12. 
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with  three  steep  peaks  which  are  clearly  related  to  the  T,  L,  and  X  valley  energies  as  shown  in 
Fig.  4.5.  The  process  of  relaxation  can  be  understood  from  the  study  of  the  time  development  of 
the  energy  distribution  which  is  shown  in  Fig.  4.6.  The  zero  of  the  time  coordinate  is  set  when 
electrons  start  to  be  injected  into  layer  1.  However,  the  actual  time  when  electrons  start  to 
enter  layer  2  is  about  2  psec  later  since  they  have  to  travel  through  layer  1  first.  As  shown  in 
the  figure,  the  electron  distribution  is  initially  very  hot  (t=2.5  psec)  while  a  peak  begins  to 
develop  near  the  T  valley  minimum  due  to  the  strong  impact  ionization.  As  the  electrons  experi¬ 
ence  increasingly  the  enchanced  electron-phonon  interactions,  the  distribution  cools  down  to 
lower  energy  and  the  peaks  near  L  and  X  valley  minima  start  to  form.  By  the  time  t=6  psec, 
electrons  have  reached  the  layer  2-layer  3  interface  and  exhibit  a  distribution  with  three  pro¬ 
nounced  peaks.  Finally,  the  peak  near  the  T  valley  minimum  grows  as  the  distribution  reaches 
the  steady  state. 

In  principle,  there  may  be  several  possible  reasons  for  these  peaks.  To  clarify  this  point,  a 
calculation  without  the  conduction  band  offset  has  been  performed.  In  this  case,  GaAs  is  used 
for  layers  1,  2,  and  3.  The  electric  field  of  300  kV/cm  is  used  for  layer  1  and  100  kV/cm  for 
layers  2  and  3.  The  results  are  plotted  in  Fig.  4.7.  Again,  the  electron  distribution  in  layer  1  is 
obtained  within  40  A  from  the  interface  while  the  distribution  in  layer  2  is  calculated  as  an 
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average  of  the  whole  layer  of  400  A  width.  The  results  for  this  case  show  that  an  overheated, 
smooth  distribution  enters  layer  2  (this  is  true  even  in  the  presence  of  a  discontinuity  as  for  Fig. 
4.5)  and  relaxes  to  a  cool,  smooth  distribution.  Steep  peaks  near  the  valley  minima  as  observed 
in  Fig.  4.5  do  not  appear.  Therefore,  the  peaks  near  the  T  and  L  valley  minima  are  due  tc  the 
effects  of  the  well.  Clearly  the  enhanced  probability  of  impact  ionization  and  the  associated 
energy  loss,  and  the  confinement  effect  due  to  the  well  boundary  (the  layer  2-layer  3  potential 
barrier)  will  lead  to  an  accumulation  of  cooler  electrons  at  the  F  valley  minimum  Indeed,  a  pro¬ 
nounced  peak  can  be  found  at  low  energies  as  shown  in  Fig.  4  5.  The  growth  of  the  peak  near 
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Figure  4.5  Electron  energy  distributions  in  structure  (b)  of  Fig.  4.1.  The  electric  field  in 
layer  1  is  500  kV/cm,  while  the  field  in  layers  2  and  3  is  100  kV/cm. 
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Figure  4.6 


Time  evolution  of  electron  energy  distribution  in  layer  2  for  structure  (b)  of  Fig 
4.1.  The  operating  conditions  are  those  of  Fig  4  5 
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Figure  4.7 


Electron  energy  distributions  in  layers  1  and  2.  GaAs  is  used  for  all  the  layers. 
The  electric  field  in  layer  1  is  300  kV/cm,  while  the  field  in  layers  2  and  3  is  100 
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the  T  minimum  in  Fig.  4.6  (between  t=6  psec  and  t=l2  psec)  is  due  to  the  further  confinement 
of  cold  electrons.  Due  to  the  potential  barrier,  the  elections  in  the  L  valley  minima  also  experi¬ 
ence  confinement.  However,  the  band  offset  for  the  L  valleys  is  very  small  (0  07  eV  when  the  A1 
mole  fraction  is  0.36),  and  the  confinement  by  this  small  barrier  alone  can  not  explain  the  dis¬ 
tinctive  peak  exhibited  in  Fig.  4.5.  These  peaks  at  higher  energies  close  to  the  L  and  X  minima 
can  be  explained  by  the  following  two  effects.  When  electrons  enter  the  quantum  well  and 
abruptly  gain  energy  from  the  band  offset,  the  distribution  experiences  the  strong  enhancement 
of  intervalley  scattering.99  Hence,  due  to  this  enhancement  and  the  differences  in  the  band 
structures,  the  L  and  X  valleys  are  overpopulated  at  the  beginning  stage.  However,  the  scatter¬ 
ing  rates  of  electrons  from  the  L  and  X  valleys  back  to  the  F  valley  are  small  due  to  the  small 
density  of  states  at  the  T  minimum  (compared  to  L  and  X).  Combined  with  the  rapid  cooling  of 
the  distribution  as  it  propagated  in  the  GaAs  layer,  the  electrons  at  the  L  and  X  valley  minima 
encounter  a  bottleneck  for  their  return  to  lower  energies.  Therefore,  as  the  overheated  electron 
gas  cools  down  fast,  very  pronounced  peaks  develop  at  the  L  and  X  minima  in  contrast  to  the 
case  of  homogeneous  structure  which  shows  only  relatively  broad  peaks  encompassing  the  ener¬ 
gies  close  to  both  L  and  X  minima. 

These  effects  occur,  of  course,  only  if  the  well  is  not  too  wide.  As  the  well  width  is 
increased,  the  energy  distribution  gradually  approaches  the  smooth  form  of  the  steady-state  dis¬ 
tribution  as  shown  in  Fig.  4.8.  For  the  case  of  the  larger  well  width,  the  cold  electrons  can  be 
reheated  by  the  relatively  high  electric  field  (100  kV/cm)  while  the  total  number  of  ionization 
events  is  virtually  constant.  Therefore,  the  population  in  the  low-energy  peak  decreases 
significantly  and  the  bottleneck  effects  diminish  also  due  to  the  longer  relaxation.  The  distribu¬ 
tion  shows  similar  changes  as  the  electric  field  applied  becomes  smaller.  Figure  4  9  shows  the 
energy  distribution  in  the  GaAs  layer  with  the  electric  field  in  layer  1  ranging  from  500  kV  cm 
to  100  kV/cm.  Although  the  distribution  entering  the  GaAs  layer  is  very  different  in  each  case, 
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electrons  in  the  GaAs  layer  have  very  similar  distribution  through  the  rapid  relaxation  process. 
Yet  the  distributions  are  not  fully  relaxed  and  show  the  dependence  on  the  initial  condition. 
The  electrons  injected  with  smaller  energy  (i.e.,  smaller  field)  have  a  less  populated  high-energy 
tail  and,  therefore,  a  smaller  T  valley  peak.  Also  the  peak  at  the  L  valley  minima  becomes  less 
pronounced  with  smaller  injection  energy  since  the  less  overheated  distribution  can  relax  faster. 
As  the  band  offset  is  decreased,  the  confinement  and  the  bottleneck  effects  are  expected  to 
decrease  rapidly.  When  the  offset  is  reduced  to  0,10  eV,  the  peaks  at  the  T  and  L  minima  virtu¬ 
ally  disappear  as  shown  in  Fig.  4.10.  It  is  also  interesting  to  note  that  the  distribution  for  the 
smaller  potential  energy  gain  (i.e.,  smaller  band  discontinuity)  has  a  larger  high-energy  tail. 
This  is  because  with  larger  band  discontinuities  the  distribution  cools  down  mo-e  due  to  the 
stronger  confinement  of  cold  electrons. 

The  above  effects  of  electron  reflection,  enhanced  intervalley  scattering  and  impact  ioniza¬ 
tion  due  to  overheating,  and  the  backscattering  bottleneck  correspond  to  the  band  structure  of 
the  material  in  a  significant  way.  Figure  4.11  shows  a  two-dimensional  plot  ((001)  plane)  of  the 
momentum  distribution  function  in  structure  (b)  of  Fig.  4.1.  In  this  figure,  the  momentum  dis¬ 
tribution  clearly  reveals  the  cubic  symmetry  of  the  crystal.  A  pronounced  accumulation  of  elec¬ 
trons  in  the  T-sphere  and  L-  and  X-  ellipsoids  makes  the  band  structure  highly  visible.  It  is 
interesting  to  note  that  the  shift  due  to  the  electron  drift  (in  the  <100>  direction)  is  rather 
small  in  spite  of  the  still  high  electric  field  of  100  kV/cm.  Therefore,  the  standard  approxima¬ 
tion  of  small  drift  in  hot  electron  theory  is  justified  However,  the  picture  of  a  displaced  Fermi- 
sphere  is  totally  incorrect.  The  constant  energy  surface  is  a  very  special  one  in  this  case,  and 
one  should  rather  call  it  a  displaced  Boltzmann  cube  The  cubic  symmetry  of  the  crystal  is  also 
found  in  the  distribution  function  for  homogeneous  GaAs.  However,  due  to  the  lack  of  accumu¬ 
lated  electrons  at  the  valley  minima,  the  structure  is  less  pronouned  and  the  distribution  of  elec¬ 
trons  is  more  diffuse,  as  shown  in  Fig.  4.12.  When  a  high  electric  field  is  applied,  the  cubic 
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Figure  4  iO  Electron  energy  distributions  m  layer  2  for  structure  (bl  of  F'lg.  41  The  electric 
field  in  layer  1  is  500  kV 'crn,  while  the  held  in  layers  2  and  3  is  100  kV  cm  The 
A1  mole  fraction  varies  from  0.36  (0.31  eV)  to  0  12  (0  10  eV) 
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Figure  4.11  Electron  momentum  distribution  in  layer  2  plotted  on  the  (001)  plane.  The 
operating  conditions  are  those  of  the  300  kV/cm  case  in  Fig.  4.9. 
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symmetry  is  not  clearly  reflected  in  the  momentum  distribution  as  shown  in  Fig.  4.13.  Since  the 
distribution  is  so  widely  scattered,  the  details  of  the  band  structure  (i.e.,  T-sphere,  L-  and  X- 
ellipsoids)  are  blurred.  The  picture  of  a  displaced  sphere  can  be  more  successfully  used  in  this 
case.  A  larger  shift  in  the  distribution  is  also  noted  in  this  case. 

We  now  turn  to  macroscopic  transport  effects,  in  particular,  to  impact  ionization.  Impact 
ionization  in  heterostructures  is  a  process  of  high  complexity  due  to  its  nonlocal  nature  and 
heavy  dependence  on  (subtle)  details  of  the  band  structure.  The  electron  impact  ionization  rate 
ot  in  the  GaAs  layer  is  studied  below  as  a  function  of  distance.  The  zero  of  the  coordinate  sys¬ 
tem  is  chosen  at  the  layer  1-layer  2  interface.  Figure  4.14  shows  the  rate  a  in  structure  (a) 
shown  in  Fig.  4.1.  In  this  figure,  o  exhibits  what  one  could  term  an  overshoot  at  a  certain  dis¬ 
tance  after  the  electrons  enter  the  GaAs  layer.  To  clarify  the  reasons  for  this  effect,  the  evolu¬ 
tion  of  the  energy  distribution  function  has  been  studied  as  the  electrons  leave  the  AlGaAs  layer 
and  transfer  to  the  GaAs  layer  experiencing  the  same  field  everywhere  (i.e.,  in  structure  (a)). 
This  study  shows  that  there  are  two  pronounced  effects  which  influence  a.  The  first  is  the  well- 
known  enhancement  of  the  ionization  due  to  the  band  edge  discontinuity  which  overheats  the 
electron  gas  and  increases,  therefore,  the  ionization  probability.12  This  effect  alone  would  lead  to 
a  monotonical  decrease  in  a  from  its  overheated  value  to  the  steady-state  value. 

The  second  effect  arises  from  the  difference  in  the  band  structure  of  the  neighboring  materi¬ 
als.  Depending  on  the  the  density  of  states  at  high  energies,  the  high-energy  tail  of  the  distribu¬ 
tion  can  be  larger  or  smaller  in  different  materials.  If  the  density  of  the  states  is  large,  the 
scattering  by  phonons  increases.  As  a  result,  the  energy  loss  becomes  high  and  the  high-energy 
tail  of  the  distribution  is  suppressed.  In  our  particular  band  structure  models,  the  density  of 
states  at  and  above  2  eV  is  larger  in  AlGaAs.  Therefore,  the  high-energy  tail  in  the  distribution 
is  smaller  in  the  AlGaAs  layer  than  in  the  GaAs  layer.  As  electrons  propagate  into  layer  2,  they 
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Figure.  4,13 


Electron  momentum  distribution  in  layer  2 
operating  conditions  are  those  of  Fig  4.2. 


plotted  on  the  (001)  plane. 
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Figure  4.14  Impact  ionization  rates  in  layer  2  for  structure  (a)  of  Fig.  4.1.  Uniform  electric 
fields  ranging  from  500  kV/cm  to  300  kV/cm  are  applied  throughout  the  layers 
(layers  1,  2,  and  3).  The  zero  of  the  coordinate  system  is  chosen  at  the  layer  1- 
laver  '2  interface. 
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abruptly  gain  energy  from  the  band  edge  discontinuity.  This  enhances  the  ionization  rate. 
However,  in  the  first  portion  of  the  well,  the  distribution  still  has  a  rapidly  decreasing  high- 
energy  tail,  which  exhibits  the  shape  typical  for  AlGaAs  as  shown  in  Fig.  4.15.  Since  the 
scattering  rate  at  high  energies  is  smaller  in  GaAs  (small  density  of  state),  electrons  can  be 
accelerated  toward  the  more  accessible  high-energy  tail  in  GaAs  fast  enough  so  that  the  ioniza¬ 
tion  rate  can  increase  while  the  overheating  of  the  distribution  is  relaxed  on  average.  Therefore, 
a  increases  over  a  distance  which  is  needed  for  the  high-energy  distribution  to  adjust  itself  to 
the  GaAs  band  structure.  Subsequently  a  decrease  to  the  steady-state  value  follows.  As  the 
electric  field  in  layer  2  is  decreased,  the  two  effects  competing  with  each  other  may  still  result  in 
a  monotonically  decreasing  a.  Even  though  electrons  begin  to  adapt  to  the  band  structure  with 
more  accessible  high  energy  states,  the  distribution  can  be  cooled  down  fast  enough  by  the  relax¬ 
ation  to  overcompensate  for  the  increase  of  the  tail.  The  ioniration  rate  for  structure  (b)  of  Fig. 
4.1  is  shown  in  Fig.  4.16.  Due  to  the  smallness  of  the  field  in  the  GaAs  layer  (100  kV/cm),  the 
relaxation  becomes  dominant  and  the  ionization  rate  monotonically  decreases  as  electrons  pro¬ 
pagate  in  this  case. 

As  mentioned  previously,  impact  ionization  is  very  sensitive  to  the  details  of  the  band 
structure.  However,  the  current  knowledge  of  the  band  structure  (especially  AlGaAs)  is  not  pre¬ 
cise.  As  a  consequence,  we  do  not  know  if  the  overshoot  of  a  indeed  exists  in  this  material  sys¬ 
tem.  However,  the  physical  reasons  for  the  effects  are  clear  and  therefore  there  may  exist  a 
material  system  which  will  show  an  effect  of  this  kind. 

Finally,  even  though  the  electron  ionization  rate  in  the  narrow  gap  layer  (GaAs)  has  been 
enhanced  by  the  heterojunction  as  shown  in  Fig  4  14,  this  enhancement  may  not  be  large 
enough  to  compensate  for  the  small  ionization  rate  in  the  wide  gap  layer  (AlGaAs)  where  virtu¬ 
ally  no  ionization  occurs.  In  spite  of  this  point,  the  a/0  ratio  can  be  still  enhanced  significantly 
The  enhancement  of  a  improves  as  the  field  becomes  low. 
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4.4  Conclusion 


A  numerical  study  of  nonlinear  transport  across  heterojunction  interface  in  strong  electric 
fields  has  been  performed.  Two  typical  heterostructures  have  been  chosen  to  analyze  the  effects 
of  a  quantum  well  on  the  electron  energy  and  momentum  distributions.  In  the  quantum-well 
APD  structure  with  a  high  field  in  the  well  layer,  electrons  experience  negligible  confinement  in 
the  well  layer  in  spite  of  heavy  ionization  and  a  potential  barrier.  On  the  other  hand,  the  doped 
quantum-well  APD  structure  (relatively  small  field  in  the  well  layer)  shows  much  stronger  relax¬ 
ation  and  confinement.  Reflections  at  the  heterointerface,  carrier  overheating  at  the  potential 
step,  and  the  L  — *•  T,  X  — *  T  back  scattering  bottlenecks  produce  pronouned  structures  in  the 
energy  and  momentum  distributions  which  directly  reveal  the  band  structure.  It  is  also  found 
that  the  ionization  rates  are  heavily  dependent  on  the  band  structure  of  the  neighboring  layer. 
An  ionization  overshoot  effect  may  exist  in  these  structures  due  to  the  band  structure  difference 
between  two  adjacent  layers. 
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CHAPTER  5 

ELECTRON  HEATING  BY  HOT  PHONON  EFFECTS  IN  HETEROSTRUCTURES 

5.1  Introduction 

Recently  the  effects  of  the  nonthermal  phonon  distributions  known  as  hot  phonon  effects 
have  received  much  attention.  The  slow  energy  relaxation  rate  of  photo-excited  hot  carriers  in 
time-resolved  measurements17,100,101  has  often  been  attributed  to  a  nonequilibrium  phonon  distri¬ 
bution  in  bulk  GaAs102  and  quantum  wells.101,103,104  Also,  the  phonon  side  band  operation  in 
photo-excited  quantum-well  heterostructure  lasers  has  been  explained  by  nonequilibrium  phonon 
effects.18,105-107  The  phonon  distribution  can  be  perturbed  as  well  when  a  large  number  of  car¬ 
riers  propagate  across  an  interface,  experience  an  abrupt  energy  gain,  and  subsequently  relax 
through  strong  phonon  emission.  As  a  general  rule,  the  phonon  distribution  will  be  disturbed 
whenever  a  significant  number  of  phonons  (especially  zone-center  phonons)  are  emitted  in  a  time 
scale  shorter  than  the  phonon  decay  time  constant. 

To  understand  the  physical  mechanism  of  this  process,  Lugli  and  coworkers  have  applied  a 
numerical  approach  and  studied  the  evolution  of  the  phonon  distribution  and  its  effects  on  the 
energy  relaxation  rate  when  hot  electrons  are  injected  by  a  strong  light  pulse.108,109  No  such 
numerical  effort  has  yet  been  made  to  analyze  the  phonon  distribution  and  its  effect  on  the  car¬ 
rier  transport  in  a  heterostructure  with  charge  carriers  propagating  across  an  interface  This 
situation  is  very  important  in  many  problems  involving  thermionic  emission  such  as  quantum- 
well  laser  diodes. 

This  chapter  presents  a  numerical  study  of  the  effects  of  nonthermal  phonon  distributions 
on  electron  transport  across  heterolayers.  For  the  microscopic  analysis,  not  only  the.  electron 
dynamics  but  also  the  phonon  dynamics  are  studied  in  detail.  In  Section  5  2  the  model  applied 
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in  the  analysis  is  briefly  described.  The  effects  of  conduction  band  discontinuity  on  ihe  phonon 
distributions  are  explained  and  thoroughly  discussed  along  with  the  influence  of  a  wide  range  of  • 

operating  conditions  in  Section  5.3.  The  conclusion  will  be  followed  at  the  end. 


5.2  Method  of  Analysis 

The  specific  structure  under  consideration  is  a  GaAs-AlGaAs  single-well  heterostructure  as 
shown  in  the  inset  of  Fig.  5.1,  which  is  essentially  the  same  with  the  structure  discussed  in 
Chapter  4.  Each  layer  is  assumed  to  be  undoped  with  an  unintentional  impurity  concentration 
of  lX1014cm-3.  The  A1  mole  fraction  used  for  layers  1  and  3  is  chosen  to  be  0.2  which  results  in 
the  conduction  band  discontinuity  of  0.17  eV  at  the  interface  assuming  AEc=0.70AEg.97 
Throughout  the  simulation,  it  has  been  assumed  that  ’ectrons  are  injected  into  the  well  layer 
after  traveling  through  layer  1.  As  discussed  in  the  previous  chapter,  layer  1  is  chosen  long 
enough  (1500  A)  to  rule  out  possible  nonlocal  effects  of  the  electron  energy  distribution.  A  400 
A -thick  GaAs  layer  forms  the  well  layer.  In  the  well  layer,  size  quantization  and  resonance 
effects  are  not  included  since  the  width  is  longer  than  the  inelastic  scattering  mean-free  path. 
For  moderately  high  electric  fields  (~  50  kV/cm),  the  scattering  mean-free  path  is  approximately 
50-100  A  63  The  specific  width  for  layer  3  is  not  important  for  this  study.  For  the  simulation, 
an  ensemble  Monte  Carlo  method  is  used  along  with  a  realistic  band  structure  in  the  GaAs 
layer.  Included  scattering  mechanisms  are  polar  optical  phonon  (POP)  and  intervalley  phonon 
scattering  as  well  as  electron-electron  scattering.  Details  of  the  treatments  (including  the  boun¬ 
dary  conditions  at  the  interface)  and  the  material  parameters  can  be  found  in  Chapter  4. 

To  account  for  the  hot  phonon  effects  the  treatment  of  POP  scattering,  which  involves  the 
zone-center  longitudinal-optical  (LO)  phonons,  has  been  modified  to  include  a  nonequilibrium 
phonon  distribution.  Although  intervalley  scattering  also  occurs  mostly  via  phonon  emission,  its 
contribution  to  the  hot  phonon  effects  is  negligible  since  intervallev  scattering  mainly  involves 
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Figure  5.1  Nonequilibrium  LO-phonon  distributions  aq(x,t)  at  three  different  locations  id  the 
well  with  the  injection  current  density  of  lX104A/cm"  The  inset  shows  the 
structure  considered  in  the  calculation  (with  a  uniform  electric  field).  The  equili¬ 
brium  phonon  distribution  n°  =  l/(exp(fTu^0/kT)-l!  is  approximately  0.34  at 
300  K.  The  ucycclization  constant  2?r/a  is  1  llX108cm  "l  for  GaAs.  The  spatial 
mesh  size  is  50  A . 
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near  zone-boundary  phonons  with  large  phonon  wavevectors.  When  the  angular  dependence  of 
the  phonon  distribution  in  momentum  space  is  neglected,  the  perturbation  in  the  phonon  dis¬ 
tribution  Anq  can  be  obtained  by  dividing  the  net  number  of  emitted  phonons  ANq  by  the  pho¬ 
non  density  of  states  Dq  (i.e. ,  Anq=ANq/Dq).  The  density  of  states  Dq  is  proportional  to  q2. 
Thus,  the  perturbation  Anq  by  phonons  with  large  q  is  suppressed  rapidly  as  q  increases  and 
POP  scattering  (which  involves  zone-center  phonons)  can  be  more  effective  in  perturbing  the 
phonon  distribution  than  intervalley  scattering  even  though  intervalley  scattering  has  a  larger 
scattering  rate.  Therefore,  only  the  perturbation  in  the  LO-phonon  distribution  by  POP  scatter¬ 
ing  is  considered.  Dving  the  simulation,  the  phonon  absorption  and  emission  events  by  POP 
scattering  are  monitored  both  in  time  and  space.  Based  upon  this  information  from  POP 
scattering  and  the  phonon-phonon  relaxation  processes  which  include  the  decay  to  zone- 
boundary  acoustic  phonons,  the  position-dependent  LO-phonon  distribution  nq(x,t)  is  recalcu¬ 
lated  every  100  fsec  as  a  function  of  the  amplitude  q  of  the  phonon  wavevector.  The  anisotropy 
of  nq(x,t)  which  may  occur  in  the  presence  of  strong  electric  fields  is  not  included  in  this  study. 
For  the  phonon-phonon  relaxation  process,  the  phenomenological  relaxation  time  rph_ph  has  been 
measured  to  be  7  psec  at  77  K,110  and  this  value  is  used  in  absence  of  direct  measurement  data 
for  300  K.  Since  the  relaxation  time  is  not  a  strong  function  of  temperature,111  rph_pl)  at 

300  K  should  indeed  be  fairly  close  to  7  psec.112  With  the  newly  calculated  phonon  distribution 
nq(x,t),  the  POP  scattering  rates  are  also  recalculated  and  the  probability  distribution  of  the 
final  states  after  POP  scattering  is  modified  to  reflect  this  new  phonon  distribution.  Possible 
effects  of  screening  have  not  been  included. 

5.3  Results  and  Discussion 

As  the  electrons  enter  the  well  from  the  barrier  (layer  1),  electrons  experience  an  abrupt 
energy  gain  at  the  heteroir forfnce.  The  overheated  electrons  then  relax  rapidly  and  cause  the 
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nonthermal  phonon  population  in  the  well  to  increase.  When  electrons  are  injected  with  higher 
energy,  more  phonons  are  generated.  However,  as  the  average  injection  energy  increases,  inter¬ 
valley  scattering  becomes  more  important  and,  thus,  many  of  the  phonons  involved  in  the  pro¬ 
cess  become  zone-boundary  phonons  which  are,  as  mentioned,  ineffective  in  disturbing  the  pho¬ 
non  distribution.  Therefore,  the  injection  of  electrons  with  average  injection  energy  well  above 
the  T-L  separation  is  not  desirable  for  the  study  of  hot  phonon  effects.  Hot  phonon  effects  are 
also  diminished  in  the  presence  of  strong  electric  fields  within  the  well  since  electrons  can  be 
quickly  reheated  and  populate  the  energy  range  of  predominant  intervalley  scattering.  In  Fig. 
5.1,  simulation  results  are  shown  for  a  uniform  electric  field  of  50  kV/cm  throughout  the  layers 
(1-3)  and  an  injection  current  density  of  lX104A/cm2  The  mesh  size  used  for  the  calculation  of 

o 

the  position-dependent  quantities  has  been  50  A  and  the  ^ero  of  the  coordinate  system  x=0  is 
at  the  layer  1-layer  2  interface  as  shown  in  the  inset  of  Fig.  5.1.  The  steady-state  LO-phonon 
distribution  is  shown  in  Fig.  5.1  at  three  different  positions  in  the  well.  The  deviation  from  the 
equilibrium  phonon  distribution  is  relatively  small  at  the  layer  1-layer  2  interface  and  grows  as 
the  electrons  travel  further  into  the  well.  This  is  because  at  the  entrance  in  the  well,  the  POP 
emission  rate  is  relatively  small  due  to  the  high  average  injection  energy  and  the  existence  of  a 
short  "near  ballistic"  range  for  the  injected  electrons.  Near  the  layer  2-layer  3  interface,  the 
reflections  due  to  the  potential  barrier  lead  to  an  accumulation  of  low-energy  electrons.  These 
electrons  with  low  energies  tend  to  emit  phonons  within  a  narrow  range  of  relatively  large 
wavevectors  as  can  be  seen  from  simple  considerations  of  energy  and  momentum  conservation.28 
Therefore,  they  cause  a  long  tail  in  the  disturbed  phonon  distribution  as  is  shown  in  the  figure. 
Similar  elongation  in  the  disturbance  of  the  phonon  distribution  has  also  been  observed  in  Ref. 
108.  These  effects  do  not  exist  in  very  wide  wells  since  the  electrons  gradually  relax  to  equili¬ 
brium  toward  the  layer  2-layer  3  interface. 
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The  nonequilibrium  phonon  distribution  affects  the  POP  scattering  rates  considerably.  The 
equilibrium  and  nonequilibrium  POP  scattering  rates  are  compared  in  Fig.  5.2.  The  rates  with 
hot  phonon  effects  in  this  figure  correspond  to  the  LO-phonon  distribution  at  the  layer  2-layer  3 
interface  in  Fig.  5.1,  where  the  phonon  perturbation  is  pronounced.  It  is  important  to  note  that 
stronger  enchancement  is  observed  in  the  absorption  rate  compared  to  the  emission  rate,  espe¬ 
cially  at  low  energies.  This  change  in  the  scattering  rate  results  in  the  difference  in  the  electron 
energy  distribution.  Figure  5.3  shows  the  electron  energy  distribution  at  the  layer  2-layer  3 
interface.  Clearly,  inclusion  of  hot  phonons  leads  to  an  increase  in  the  electron  temperature  as 
shown  by  the  fact  that  there  are  more  electrons  at  high  energies.  Due  to  the  strong  enhance¬ 
ment  in  the  phonon  absorption  rate,  cold  T  valley  electrons  in  the  well  can  reabsorb  the  non¬ 
equilibrium  phonons  which  have  been  emitted  by  the  hot  electrons  injected  over  the  barrier. 
Thus,  with  the  inclusion  of  hot  phonons,  more  electrons  transfer  to  the  L  valleys  and  contribute 
to  a  more  pronounced  peak  at  the  L  minima.  It  is  interesting  to  note  that  the  backscattering 
bottleneck  effect113  is  also  observed  at  the  energy  of  the  L  valley  minima.  Another  interesting 
point  is  that  two  different  procedures  which  have  been  used  to  calculate  the  hot  phonon  effects 
give  almost  identical  energy  distributions.  Procedure  (A)  is  the  method  described  above  in 
which  both  the  POP  scattering  rate  and  the  probability  distribution  of  the  final  states  after 
POP  scattering  is  updated.  In  procedure  (B),  the  POP  scattering  rate  is  recalculated  based  on 
the  updated  phonon  distribution,  but  the  probability  distribution  of  the  final  states  for  thermal 
equilibrium  is  used  to  choose  the  final  states  after  POP  scattering.  Evidently,  procedure  (B)  has 
a  different  preference  for  the  phonon  wavevectors  involved  in  the  scattering  and  generally  results 
in  different  final  electron  states  compared  to  procedure  (A).  However,  since  a  constant  phonon 
energy  is  used  for  POP  scattering,  the  change  in  electron  energy  is  the  same  in  both  procedures. 

The  heating  by  the  reabsorption  of  hot  phonons  reduces  the  confinement  of  electrons  in  the 
well.  Figure  5.4  shows  the  time  evolution  of  the  number  of  electrons  in  the  well  (normalized  to 


Calculated  POP  scattering  rates  in  equilibrium  and  nonequilibrium  cases  The 
nonequilibrium  rates  correspond  to  the  LO-phonoD  distribution  at  the  layer  2- 
layer  3  interface  in  Fig.  5  1. 


Energy  Distribution  (arbitrary  units) 


Energy  (eV) 


Figure  5.3  Electron  energy  distributions  at  the  layer  2-layer  3  interface.  For  procedure  (A) 
the  probability  distribution  for  the  final  states  is  updated  to  reflect  the  nonequili¬ 
brium  distribution,  while  the  probability  distribution  at  thermal  equilibrium  is 
used  for  procedure  (B). 


Number  of  Electrons  in  the  Well  (arbitrary  units) 


Figure  5.4  Time  evolution  of  the  number  of  electrons  in  the  well  (normalized  to  the  steady 
state  value)  after  the  injection  current  is  turned  off  at  t=0  The  inset  shows  the 
corresponding  development  of  the  averse  energy  in  the  well.  The  number  of 
electrons  and  the  average  energy  are  averaged,  over  layer  2. 
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the  steady-state  value)  after  the  injection  current  is  turned  off  at  t=0.  Due  to  the  residual  elec¬ 
trons  in  layer  1,  the  response  to  the  switch-off  is  delayed.  As  can  be  seen  in  the  figure,  electrons 
are  emitted  from  the  well  about  two  times  faster  than  without  the  inclusion  of  hot  phonon 
effects.  The  corresponding  development  of  the  average  electron  energy  in  the  well  with  time  is 
shown  in  the  inset  of  Fig.  5.4.  Due  to  the  emission  of  electrons  with  high  energies  above  the 
barrier,  the  average  energy  falls  rapidly  and  then  approaches  equilibrium  in  the  well.  As 
expected,  the  difference  in  average  energies  (with  and  without  the  inclusion  of  hot  phonon 
effects)  decreases  with  time  but  is  still  sizable  after  several  picoseconds  due  to  the  presence  of  hot 
phonon  absorption. 

As  the  electric  field  in  the  well  decreases,  the  electrons  experience  stronger  relaxation 
toward  lower  energies  which  results  in  an  energy  distribution  close  to  equilibrium.  The  decrease 
in  the  average  injection  energy  also  leads  to  a  reduction  in  the  electron  temperature  in  the  well. 
As  the  electron  distribution  cools  down,  the  effect  of  phonon  reabsorption  from  the  nonequili¬ 
brium  distribution  grows  since  an  increasing  fraction  of  electrons  resides  near  or  below  the 
threshold  energy  for  POP  emission  where  only  absorption  is  allowed.  Figure  5.5  shows  the  elec¬ 
tron  energy  distribution  in  the  well  when  the  electric  field  is  zero  in  all  layers,  which  is  analogous 
to  the  operating  conditions  in  quantum-well  laser  diodes.  In  contrast  to  Fig.  5.4,  the  injection  of 
electrons  into  the  well  layer  is  started  at  t=0  from  the  layer  1-layer  2  interface  by  thermionic 
emission  in  this  calculation.  The  same  injection  current  density  as  for  Fig.  5.1  is  assumed  Since 
the  simulation  is  very  time-consuming  under  these  circumstances,  electrons  are  simulated  only 
for  the  first  12  psec.  The  distribution  displayed  in  Fig.  5.4  represents  an  average  over  the  whole 
layer  2  (400A)  at  t= 12  psec.  As  shown  in  the  figure,  a  significant  portion  of  the  electrons  near 
the  threshold  energy  for  emission  absorbs  phonons  and  moves  to  higher  energy.  About  30% 
increase  in  the  electron  lein^nature  is  therefore  observed  (inset  of  Fig  5.5)  compared  to  about 
15%  increase  in  the  previous  case  (inset  of  Fig.  5.4).  The  hot  phonon  distribution  will  have  an 
even  larger  effect  in  steady  state  since  the  effect  is  still  growing  with  time  at  t= 1 2  psec. 


Energy  Dis rr ibution  (arbitra 


Figure  a 


Energy  (eV) 


o  Electron  energy  distribution  in  the  well  at  t=  12  psec  with  no  electric  field 
applied  The  time  development  in  the  average  energy  id  the  well  is  shown  in  the 
inset.  The  injection  nMrrent  is  turned  on  at  t— n  The  distribution  and  the  aver 
age  energy  are  averaged  over  layer  2. 
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The  injection  current  density  is  an  important  parameter  as  large  current  densities  enhance 
hot  phonon  effects.  A  smaller  but  still  considerable  disturbance  of  the  phonon  distribution  is 
seen  in  Fig.  5.6  with  an  injection  current  density  of  lX103A/cm?.  With  very  high  injection 
currents,  the  perturbation  of  near  zone-boundary  phonons  may  also  become  important.  The  hot 
phonon  effects  will,  of  course,  be  more  visible  at  lower  temperature.  Due  to  the  small  phonon 
occupation  at  thermal  equilibrium,  phonon  absorption  is  possible  only  in  nonequilibrium  situa¬ 
tions  at  low  temperature. 

5.4  Conclusion 

The  effects  of  the  nonequilibrium  phonon  distribution  on  electron  transport  across  heteroin¬ 
terfaces  have  been  studied  by  an  ensemble  Monte  Carlo  simulation  method.  A  strongly  dis¬ 
turbed  phonon  distribution  is  observed  as  a  consequence  of  the  conduction  band  discontinuity 
The  hot  phonon  distribution  causes  in  turn  more  energetic  electron  distribution  because  of  pho¬ 
non  reabsorption  and  reduces,  therefore,  the  confinement  of  electrons  in  the  well.  As  the  average 
electron  energy  decreases  in  the  well,  the  effect  of  phonon  reabsorption  gains  significance  since 
more  electrons  are  near  or  below  the  threshold  energy  for  POP  emission  where  only  absorption  is 


allowed. 
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CHAPTER  6 

■ 

ELECTRON  TRANSPORT 

IN  TUNNELING  HOT  ELECTRON  TRANSFER  AMPLIFIERS 


■ 


t 


6.1  Introduction 

The  idea  of  ballistic  transport  in  thin  semiconductor  layers  has  received  considerable  atten¬ 
tion.  However,  experimental  verifications  of  such  transport  have  been  rare  and  indirect. 
Recently,  Levi  et  al.15  and  Heiblum  et  al.16  have  independently  reported  direct  observation  of 
ballistic  transport  with  the  use  of  hot  electron  transistors  as  an  electron  energy  spectrometer. 
Based  on  studies  of  the  collector  current-voltage  characteristics,  they  have  deduced  the  energy 
spectra  of  electrons  which  travel  through  the  heavily  doped  short  base  region,  and  have  shown 
the  existence  of  a  ballistic  (or  quasiballistic)  peak.  Since  then,  a  number  of  theoretical  ana¬ 
lyses114-117  have  been  performed  to  investigate  the  details  of  this  nonequilibrium  transport 
phenomenon  and  to  assess  the  feasibility  of  hot  electron  spectroscopy.  Most  of  these  efforts  have 
been  devoted  to  studying  transport  in  the  base  as  a  function  of  the  collector  barrier  height.  For 
the  spectrometer  with  a  planar  doped  barrier  structure,15  such  analyses  are  appropriate  since  the 
collector  current  is  controlled  mainly  by  the  collector  barrier  height  (collector-base  bias).  For  the 
tunneling  hot  electron  transfer  amplifier  (THETA)  structure,16  however,  the  collector  current  can 
in  principle  be  influenced  not  only  by  the  collector-base  bias  but  also  by  the  details  of  transport 
in  the  relatively  thick  (~  1000  A  )  collector  barrier  region  with  a  retarding  field. 

In  this  chapter,  the  electron  transport  in  a  THETA  structure  is  studied  numerically.  The 
evolution  of  the  electron  energy  distribution  (width  and  location  of  the  peak)  in  the  base  and  the 
collector  barrier  is  simulated  in  detail  to  examine  the  validity  of  hot  electron  spectroscopy.  The 
device  characteristics  as  a  transistor  including  transfer  ratio  and  transit  time  are  also  investi¬ 
gated.  The  simulation  results  are  compared  to  the  experiments  of  Heiblum  and  coworkers.16 
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6.2  Computational  Model 

The  conduction  band  profile  of  the  structure  under  consideration  is  schematically  shown  in 
Fig.  6.1.  The  structure  consists  of  two  undoped  AlGaAs  layers  between  three  heavily  doped 
GaAs  layers  (Np  =  IX  1018  cm-3,  where  Nq  represents  the  donor  concentration).  The  A1  compo¬ 
sition  and  the  widths  of  the  two  AlGaAs  layers  are  0.3,  120  A,  and  1000  A,  respectively.  These 
conditions  are  chosen  to  match  the  experiment  by  Heiblum  et  al.16  A  slightly  larger  base  (400 
A)  has  been  chosen  to  encompass  the  width  of  the  grading  region.  The  compositional  grading 
itself  is  not  included  in  the  simulation.  For  the  electron  transport  study,  an  ensemble  Monte 
Carlo  method  has  been  used.  Analytical  T-L-X  band  structures  with  nonparabolicity  are  incor¬ 
porated  for  both  GaAs  and  AlGaAs.  Using  AE,.=0.68AEg  for  the  conduction  band  discon¬ 
tinuity,97  the  particular  Al  mole  fraction,  0.3,  corresponds  to  a  bandedge  discontinuity  of  0.25 
eV  at  the  layer  boundaries.  Material  parameters  that  have  been  used  are  described  in  Chapter 
4 

Included  scattering  mechanisms  are  polar  optical  phonon  and  intervalley  phonon  scattering 
for  the  undoped  AlGaAs  layers  (barriers),  and  elastic  ionized-impurity  and  coupled 
plasmon/phonon  scattering  as  well  as  intervalley  phonon  scattering  for  the  heavily  doped  GaAs 
layer  (base).  Due  to  the  degenerate  electron  concentration,  the  polar  optical  phonon  modes  in 
the  base  are  strongly  coupled  to  the  collective  plasmon  mode  and  the  electron-hole  excitation 
through  long-range  polarization  fields.  The  scatterings  by  the  coupled  plasmon/phonon  modes 
are  basically  treated  as  discussed  in  Chapter  2.  In  addition  to  the  three-dimensional  treatment 
in  Chapter  2,  a  quasi-twodimensional  correction  is  also  included  in  this  study.  As  the  device 
dimension  shrinks,  the  plasmon  mode  is  strongly  influenced  by  the  reduced  dimensionality  The 
electron  gas  does  not,  of  course,  support  a  collective  motion  perpendicular  to  the  layer  interface 
whose  wavelength  is  larger  than  the  layer  width  The  exact  inclusion  of  such  effects  involves  a 
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Figure  6.1  Conduction  band  profile  for  an  AlGaAs/GaAs  tunneling  hot  electron  transfer 
amplifier  when  a  retarding  field  is  applied  in  the  collector  barrier.  The  doping 
concentration  in  the  n+-GaAs  layers  is  1X1018  cm-3  The  AlGaAs  layers  are 
undoped.  The  A1  composition  is  0.3. 


102 


great  deal  of  complexity  and  has  not  been  included  in  this  study.  Instead,  a  maximum 
wavelength  limit  is  simply  imposed  to  the  plasmon-like  mode  based  on  the  base  width  (device 
limit).  The  possibility  of  interaction  with  two-dimensional  plasmons  which  propagate  parallel  to 
the  layer  interface  still  exists.  Since,  however,  the  electrons  entering  the  base  initially  have  negli¬ 
gible  parallel  momentum,  the  coupling  to  the  two-dimensional  plasmon  can  be  neglected  due  to 
the  tight  restriction  in  the  energy  and  momentum  conservation  as  discussed  in  the  following. 
Our  treatment  poses  as  a  first-order  approximation  to  include  two-dimensional  effects. 

Let  k  be  the  momentum  of  initial  state  and  k'  the  momentum  of  final  state  after  scatter¬ 
ing.  Then  the  energy  conservation  for  plasmon  emission  gives 


ti2 

2m* 


k'  2  +  lTwp(q)  =  Ji-k 
2m 


2 


(6.1) 


when  a  parabolic  band  structure  is  assumed  with  the  effective  mass  m*.  The  plasmon  dispersion 
relation  in  two  dimensions118  gives 

w2(q)  ~  aq  +  0q2  (6.2) 

where  a  =  2e2Ep/B2  and  /?  ~  3Ep/2m*.  Here,  Ep  is  the  Fermi  energy.  Assuming  k— ’ q=k' ,  Eq. 
(6.1)  gives 


_q_ 

2k 


(6.3) 


Therefore,  for  a  given  k,  there  is  a  critical  angle  ©c  above  which  an  appropriate  q  (in  turn,  final 
state)  can  not  be  found.  Assuming  aq»/?q2  which  is  appropriate  at  long  wave  length,  an  esti¬ 
mate  of  ©c  ~  40  °  is  obtained  with  m*/m0  =  0.067  and  E(lT)  =  0.4  eV  for  an  electron  density  of 
lxl012cm~2  at  T  =  0  K  As  E(lT)  decreases,  ©,,  also  decreases  to  a  smaller  value.  This  estimate 
shows  that  two-dimensional  plasmon  scattering  is  an  ineffective  scattering  mechanism  unless  the 
electrons  have  considerable  parallel  momentum  (i.e.,  0<©<9C).  Although  the  conservation  of 


perpendicular  momentum  (TT^  =  Ic^' )  can  not  be  justified  in  full  due  to  the  lack  of  translational 
invariance,  this  estimate  shows  qualitatively  that  two-dimensional  plasmon  scattering  is  of  lim¬ 
ited  importance. 

Some  of  the  calculated  inelastic  scattering  rates  are  shown  in  Fig.  6.2  as  a  function  of 
energy.  Two  different  procedures  have  been  used  to  obtain  the  total  inelastic  scattering  rates. 
Procedure  (A)  includes  the  device  limit  (base  width)  in  the  plasmon-like  mode  as  described 
above,  while  procedure  (B)  does  not.  As  can  be  seen  in  the  figure,  these  two  procedures  exhibit 
sizable  differences  in  the  scattering  rate.  Although  the  Landau  damping  mode  is  the  most  dom¬ 
inant  scattering  mechanism  and  is  treated  independent  of  the  base  width  in  this  study,  the 
difference  between  the  two  treatments  (A  and  B  above)  of  the  plasmon-like  mode  affects  the  out¬ 
come  of  the  transport  study  to  a  considerable  degree.  This  further  emphasizes  the  importance  of 
precise  treatment  of  the  finite  base  size.  For  ionized  impurity  scattering,  the  well-known 
Brooks-Herring  model  has  been  used  along  with  Thomas-Fermi  screening.28  After  the  scattering, 
electron  states  are  randomly  chosen  based  on  energy  and  momentum  conservation.  Due  to  the 
broadened  resonance  frequencies,  the  Landau  damping  mode  involves  more  difficulties  in  the 
selection  of  final  states  than  other  mechanisms.  In  this  study,  the  amount  of  energy  loss  So;  for 
a  particular  scattering  event  by  the  Landau  damping  mode  is  determined  by  a  precalculated  dis¬ 
tribution  F(cj)  shown  in  Fig.  6.3.  After  choosing  fTw,  the  change  in  momentum  is  decided  by 
another  look-up  table.  The  final  state  of  the  electron  excited  from  the  Fermi  sea  is  subsequently 
picked  based  on  the  changes  in  the  scattered  electron.  The  Pauli  exclusion  principle  (zero  tem¬ 
perature)  is  included  for  the  final-state  selection. 

Since  the  temperature  considered  is  very  low  (4.2  K),  it  is  assumed  that  electrons  are 
injected  into  the  base  mainly  by  tunneling  through  the  emitter  barrier.  Assuming  further  that 
the  emitter  barrier  tunneling  probability  D(E)  depends  only  on  longitudinal  electron  energy  Ex, 
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Figure  6  2 


Calculated  inelastic  scattering  rates  for  T  valley  electrons  in  the  base.  Procedure 
(AJ  for  total  inelastic  scattering  rate  includes  the  device  limit  (base  width)  m  the 
plasmon-like  mode,  while  procedure  (B)  does  not. 


F(oj)  (arbitrary  units) 
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Figure  6.3  Probability  distribution  F(w)  for  the  amount  of  energy  loss  iTw  in  a  particular 
scattering  event  by  the  Landau  damping  mode.  The  distribution  F(w)  is  defined 

3  max 

as  f  Im  —  *  .  },  and  is  calculated  as  described  in  Chapter  2.  Nonparaboli- 

q L  * 

city  is  included  in  the  calculation.  To  obtain  the  probability  distribution  for  a 
specific  initial  energy  E,  F(w)  should  be  renormalized  to  the  area  under  F(w)  from 
fTw  =  0  to  E  —  Ep. 


it  can  be  shown  that  the  energy  distribution  of  the  injected  electrons  is  proportional  to 
(Ep—  EX)D(EX)  for  Ex<Ep  at  0  K,  where  Ep  is  the  Fermi  energy  in  the  emitter.  In  the  simula¬ 
tion,  electrons  are  continuously  launched  at  every  10  fsec  from  the  emitter  barrier-base  interface 
based  on  the  energy  distribution  discussed  above,  and  are  tracked  until  they  leave  the  device  or 
join  the  cold  base  electron  sea.  The  cold  base  electrons  are  not  simulated  unless  excited  from  the 
Fermi  sea  through  scatterings  by  the  Landau  damping  modes.  The  transfer  matrix  formalism 
discussed  in  Chapter  2  is  used  in  the  calculation  of  tunneling  probability  D(EX).  Parallel 
momentum  is  initially  neglected  since  it  is  much  smaller  than  the  perpendicular  momentum. 
Coherent  interference  effects39  are  not  considered  in  this  study.  At  the  interfaces,  electron 
transfer  for  the  L  and  X  valley  electrons  is  treated  classically  by  conserving  energy  and  parallel 
momentum.  An  electron  can  cross  the  interface  when  it  has  high  enough  energy  Ex  to  overcome 
the  barrier.  The  valley  to  which  the  electron  belongs  is  assumed  to  be  unchanged  during  the 
transfer.  For  the  T  valley  electrons,  quantum  mechanical  reflection  is  included.  The  reflections 
at  the  interfaces  are  assumed  to  be  specular.  When  an  electron  in  the  base  has  smaller  energy 
than  the  interface  barrier  height,  it  is  assumed  that  the  electron  either  reaches  the  base  electrode 
or  recombines  with  a  hole.  Since  the  continuous  launching  scheme  is  used  as  mentioned  above, 
the  collector  current  can  be  obtained  by  the  number  of  electrons  that  reach  the  collector  region 
per  unit  time.  Finally,  the  effect  of  band  bending  is  not  included  in  the  simulation. 

6.3  Results  and  Discussion 

The  electron  energy  distributions  are  shown  in  Fig  6.4  for  two  different  locations  and  an 
emitter-base  bias  Veb  of  -0.35  V  and  a  collector- base  bias  Vcb  of  -0.10  V  The  injected  energy 
distribution  is  estimated  as  described  above  and  shows  a  peak  at  ~  15  meV  below  the  Fermi 
level  with  a  full  width  of  ~  35  meV  at  half  maximum.  The  distribution  after  the  base  is 
obtained  at  the  collector  barrier  side  of  the  interface  based  on  the  incoming  electrons  from  the 
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Collector-Base  Bias  VCb  (Volts) 

0  -0.05  -0.  0  -0.15  -0.20 


Figure  6.4  Electron  energy  distributions  at  various  locations  along  with  3Ie/<9Vcb  when  the 
emitter-base  bias  is  -0.35  V.  The  electrons  are  ballistioallv  inie-'ted  into  the  base 
from  the  emitter  through  tunneling.  The  distribution  after  the  base  has  been 
obtained  at  the  collector  barrier  side  of  the  base-collector  barrier  interface  based 
on  the  incoming  electrons  from  the  base.  The  zero  of  the  potential  energy  is 
chosen  at  the  beginning  of  the  collector  barrier  which  is  AEC  above  the  base 
potential.  <9Ic/<9Vcb  has  been  calculated  directly  from  the  simulation  result  of  the 
collector  current-voltage  characteristics. 
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base.  The  electrons  traveling  back  to  the  base  are  not  considered  in  the  distribution.  Since  the 
electrons  which  are  backscattered  into  the  base  from  the  collector  barrier  can  reenter  the  collec¬ 
tor  barrier  after  suffering  further  scattering  in  the  base,  the  distribution  after  the  base  is  a  func¬ 
tion  of  Vcb.  This  dependence,  however,  is  not  prominent  in  this  study  due  to  the  relatively  large 
conduction  band  discontinuity  AEC.  Hence,  the  distribution  calculated  at  a  particular  Vcb  can 
reflect  the  intrinsic  base  transport  characteristics.  The  zero  of  the  potential  energy  is  chosen  at 
the  beginning  of  the  collector  barrier  which  is  AEC  above  the  base  potential.  As  can  be  seen  in 
the  figure,  the  electrons  preserve  the  shape  of  the  injected  distribution  after  traveling  through 
the  base.  While  the  width  is  slightly  widened  by  the  scattering,  the  peak  of  the  distribution 
after  the  base  has  not  been  displaced,  thus  demonstrating  the  ballistic  (or  quasiballistic)  tran¬ 
sport  in  the  base.  Indeed,  a  significant  portion  (~  55%)  of  the  electrons  which  enter  the  collec¬ 
tor  barrier  has  passed  through  the  base  without  a  collision.  Some  electrons  which  suffer  small 
angle  scatterings  (especially  elastic)  overcome  the  potential  barrier  and  contribute  to  the  long  tail 
in  the  distribution.  Contributions  by  the  electrons  excited  from  the  Fermi  sea  is  negligible.  No 
distinctive  secondary  peak  due  to  phonon  replica  is  observed  in  the  tail.  As  discussed  previ¬ 
ously,  the  scattering  by  the  Landau  damping  mode  is  the  most  frequent  mechanism  in  the  base. 
However,  the  T  — ►  L  intervalley  scattering  is  more  detrimental  to  the  base  transport  due  to  its 
randomizing  nature.88  Electrons  which  suffer  the  intervalley  scattering  in  the  base  have  little 
chance  to  overcome  the  barrier  and  eventually  contribute  to  the  base  current. 

To  investigate  the  experimental  evaluation  of  the  distribution,  we  have  calculated  the 

derivative  3Ic/3Vcb  and  plotted  it  also  in  Fig.  6.4.  It  is  seen  that  <9Ic/dVcb  exhibits  a  pro¬ 
nounced  peak  which  can  be  clearly  related  to  the  ballistic  peak  of  the  distribution  after  the  base 

region  (thus,  the  injected  distribution)  The  location  and  the  width  of  the  peaks  obtained  from 

these  three  different  curves  are  identical  within  the  simulation  error.  Such  coincidence  implies 
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that  at  sufficiently  large  negative  Vcb  only  ballistic  (or  quasiballistic)  electrons  contribute  to  the 
collector  current,  and  the  distribution  at  the  end  of  the  collector  barrier  still  resembles  closely 
the  injected  distribution.  This  is  clearly  shown  in  Fig.  6.5  where  the  energy  distribution  at  the 
collector  barrier-collector  interface  is  plotted  for  various  Vcb.  For  Vcb  =  Veb  -I-  AEC  (here,  Vcb  = 
-0.10  V),  the  distribution  reproduces  the  shape  of  the  injected  distribution  after  traveling 
through  the  base  and  the  collector  barrier.  More  than  75%  of  the  transmitted  electrons  are 
ballistic  in  th!~  caoe.  As  the  negative  bias  increases,  the  electrons  keep  virtually  the  same  distri¬ 
bution,  losing  only  the  portion  which  is  reflected  by  the  potential  barrier.  A  slight  growth  in  the 
peak  is  also  observed,  which  is  due  to  reduced  scattering  in  the  collector  barrier.  These  results 
prove  clearly  the  concept  of  ballistic  transport  and  its  detection  in  THETA  devices.  As  Vcb 
becomes  less  negative,  the  bulk  of  the  scattered  electrons  start  to  enter  the  collector  Since  the 
collector  current  in  this  regime  is  strongly  influenced  by  the  details  of  scattering  in  the  collector 
barrier,  it  is  difficult  to  deduce  information  on  the  base  transport  from  dl e/<9Vci,.  The  evolution 
of  the  energy  distribution  and  <9Ic/dVcb  have  also  been  studied  at  Veb  =  -0.40  V.  Due  to  the 
increase  in  injection  energy,  electrons  are  subject  to  stronger  scattering  in  both  the  base  and  the 
collector  barrier.  However,  no  qualitative  deviation  from  Fig.  6.4  is  observed  as  shown  in  Fig. 
6.6. 

Figure  6.7  illustrates  the  ballistic  transfer  ratio  in  the  collector  current  as  a  function  of  Vcb 
for  Veb  =  -0.35  V  and  -0.40  V.  At  large  negative  Vcb,  the  ballistic  transfer  ratio  is  large  and 
varies  slowly  with  Vcb.  Perfect  ballistic  transport  is  not  realized  even  at  very  high  bias  due  to 
the  elastic  small  angle  scattering  by  ionized  impurities  in  the  base.  A  drastic  decrease  in  the 
ballistic  ratio  follows  as  Vcb  becomes  less  negative.  Since  the  massive  contribution  of  scattered 
electrons  to  the  collector  current  causes  this  decrease,  the  location  of  the  sudden  drop  in  the 
ballistic  transfer  ratio  for  Veb  =  -0.35  V  matches  well  with  the  increase  in  <9Ic/<9Vcb  shown  in  Fig. 


■ 


Energy  Distribution  (arbitrary  units) 


Figure  6.5  Energy  distributions  at  the  collector  barrier  side  of  the  collector  barrier-collector 
interface  with  various  collector-base  bias  conditions  wheD  the  emitter-base  bias  is 
-0.35  V.  Each  curve  is  normalized  differently  so  that  the  ratio  of  two  distribu¬ 
tions  at  a  given  energy  represents  the  ratio  of  collector  current  with  that  energy 
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at  various  locations  along  with  d\J <9Vcb  when  the 
The  distribution  after  the  base  is  obtained  at  the 
Other  conditions  are  the  same  with  Fig.  6.4. 
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Figure  6  7 


Ballistic  transfer  ratio  in  the  collector  current  as  a  fuDctioD  of  collector- base  bias. 
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6.4,  as  discussed  previously.  With  an  increase  in  injection  energy,  ballistic  transport  becomes 
less  probable  for  a  given  Vcb  due  to  the  increased  intervalley  scattering  rate. 

Thus,  the  results  of  our  numerical  study  are  in  excellent  overall  agreement  with  the  experi¬ 
ment  by  Heiblum  et  al.16  However,  there  are  points  that  need  consideration  for  a  complete  com¬ 
parison.  In  the  experiment,  the  width  of  the  distribution  obtained  from  <?Ic/<9Vcb  (~  60  meV)  is 
larger  than  in  the  simulation  (~  35  meV).  Most  of  the  broadening  comes,  in  our  opinion,  from 
the  band  bending  at  the  emitter-emitter  barrier  interface.  Due  to  the  formation  of  the  accumu¬ 
lation  region,  the  injected  distribution  of  the  tunneling  current  can  be  widened  considerably.53 
Similar  band  bending  occurs  at  both  the  interfaces  of  the  base  region  and  may  affect  the  electron 
transport  considerably.  To  investigate  such  effects,  a  one-dimensional  analysis  as  described  in 
Chapter  2  can  be  used.  Due  to  the  finite  base  resistance,  the  choice  of  base  potential  should  be 
considered  carefully  when  the  base  width  is  very  small.  The  reflections  at  the  collector  contact, 
which  is  less  than  ideal,  can  also  broaden  the  distribution  in  the  experiment.  The  experimental 
result  for  5Ic/9Vcb  shows  that  as  Vcb  becomes  less  negative,  a  smaller  and  relatively  flat  distribu¬ 
tion  follows  the  ballistic  peak  when  compared  to  the  numerical  result  plotted  in  Figs.  6  4  and 
6.6.  This  suggests  that  the  simulation  predicts  a  larger  population  of  electrons  at  low  energies 
and,  thus,  a  stronger  energy  loss  in  the  direction  of  propagation  by  scattering  than  the  experi¬ 
ment  includes.  Considering  that  the  base  transfer  ratio  agrees  closely  in  both  cases,  the  scatter¬ 
ing  rate  in  the  collector  barrier  (mostly  polar  optical  phonon)  might  have  been  slightly  overes¬ 
timated.  In  the  simulation,  about  36%  of  the  injected  electrons  transfer  the  base  ballistically 

In  normal  transistor  operation,  positive  Vcb  is  applied  to  obtain  high  Ic/Ie  ratio.  Due  to 
strong  intervalley  scattering  in  the  collector  barrier,  ballistic  transport  through  the  device  is 
highly  improbable  with  positive  Vcb.  However,  electrons  still  travel  significant  distances  ballisti¬ 
cally  and,  thus,  have  a  higher  velocity  than  the  saturation  velocity  Figure  6.8  shows  the  total 
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Figure  6.3  Total  transit  time  as  a  function  of  collector- base  bias 
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transit  time  as  a  function  of  Vcb.  There  exists  a  negative  Vcb  range  which  minimizes  the  transit 
time,  corresponding  to  the  condition  of  large  ballistic  probability  and  simultaneously  relatively 
little  deceleration  of  the  electrons  by  the  applied  retarding  field.  As  Vcb  becomes  less  negative, 
the  transit  time  increases  sharply  at  the  onset  of  intervalley  scattering  and  gradually  approaches 
to  a  saturation  value.  At  Vcb  =  1.0  V,  the  transit  time  is  approximately  0.71  psec  and  0.74  psec 
for  Veb  =  -0.35  V  and  -0.40  V,  respectively,  while  the  base  transit  time  is  approximately  0.05 
psec  and  0.06  psec.  It  is  difficult  to  define  a  base  transit  time  with  negative  Vcb  due  to  the  back- 
scattering  of  electrons  from  the  collector  barrier.  Without  considering  the  backscattered  elec¬ 
trons  which  hardly  contribute  to  the  collector  current,  the  base  transit  time  for  negative  Vcb  is 
expected  to  be  close  to  the  value  with  positive  bias.  The  total  transit  time  is  the  parameter 
which  gives  the  ultimate  limit  of  the  switching  speed.  To  reduce  the  switching  time,  a  short 
base  width  is  desired.  However,  a  short  base  width  brings  large  base  resistance,  which  also  lim¬ 
its  the  switching  speed.  Higher  doping  in  the  reduced  base  may  be  the  solution  for  this  problem. 
While  keeping  the  base  resistance  low,  higher  doping  in  the  narrow  base  may  cause  considerably 
smaller  increase  in  the  scattering  rate  of  the  injected  electrons  than  in  the  bulk  due  to  the 
reduced  dimensionality  in  the  density  of  the  final  states.  Figure  6.9  shows  the  collector  current 
switching  characteristics  as  Vcb  switches  from  -1.0  V  to  1.0  V  with  Veb  =  -0.35  V.  At  t  =  0.75 
psec,  the  collector  current  reaches  to  more  than  90%  of  the  steady-state  value,  which  shows  good 
agreement  with  the  total  transit  time.  There  is  a  possibility  that  the  current  rise  time  can  be 
much  shorter  than  the  transit  time.  When  a  small  negative  Veb  is  applied  for  the  off  state, 
mobile  electrons  can  be  effectively  stored  in  the  collector  barrier  and  contribute  to  the  collector 
current  rapidly  when  the  bias  is  switched  to  positive.  An  overshoot  in  current,  is  expected  in  this 
situation.  As  the  current  density  increases,  the  band  bending  in  the  collector  barrier  due  to 
mobile  electrons  becomes  more  important.  The  effect  of  a  self-consistent  field  has  been  found  to 
be  negligible  even  at  500  A/cm2.  Notice,  however,  large  current  density  can  also  bring  the 
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Figure  8.9  Collector  current  switching  characteristics  as  the  collector- base  bias  switches  from 
-1.0  V  to  1.0  V  at  t  =  0  with  the  emitter-base  bias  of  -0.35  V.  The  collectoi 
current  is  normalized  to  the  emitter  current. 
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nonequilibrium  phonon  distribution  known  as  hot  phonon  effect101  in  the  base  and  can  alter  the 
base  transport  significantly. 

6.4  Conclusion 

I  •  a  Q 

In  this  chapter,  the  electron  transport  in  a  THETA  device  with  a  base  width  of  400  A  and 

o 

the  collector  barrier  width  of  1000  A  has  been  numerically  studied  at  4.2  K.  Using  an  ensemble 
Monte  Carlo  method  which  includes  coupled  plasmon-phonon  interaction  and  the  Pauli  exclusion 
principle  in  the  heavily  doped  base  region,  we  have  shown  the  existence  of  nearly  ballistic  trans¬ 
port  in  the  base  and  the  collector  barrier  for  a  negative  Vcb.  From  dIc/3Vcb,  the  ballistic  peak  of 
the  distribution  after  the  base  (thus,  the  injected  distribution)  has  been  recovered  which 
confirms  the  ideas  of  hot  electron  spectroscopy  in  THETA  structures  by  Heiblum  and  coworkers. 
With  positive  Vcb,  ballistic  transport  through  the  device  is  highly  improbable  due  to  strong 
|  intervalley  scattering  in  the  collector  barrier.  However,  electrons  still  travel  much  of  the  device 

ballistically  and,  thus,  have  a  higher  velocity  than  the  saturation  velocity.  The  total  transit 
time  which  gives  the  ultimate  limit  of  the  switching  speed  is  approximately  0.71  psec  at  Vcb  = 
I  -1.0  V  when  Veb  =  -0.35  V.  To  obtain  optimum  transit  time  and  current  gain,  the  electron 

injection  energy  in  the  base  should  be  lower  than  the  threshold  for  T  — ►  L  phonon  emission. 
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CHAPTER  7 


SUMMARY 

Electron  transport  in  III- V  semiconductors,  especially  the  GaAs/AlGaAs  material  system,  is 
studied  in  various  nonequilibrium  situations.  Throughout  the  study,  a  Monte  Carlo  simulation 
method  is  used  for  the  first-principle  analysis  of  transport  properties  under  the  semiclassical 
Boltzmann  transport  picture.  Physical  models  to  include  such  complicated  effects  as  coupled 
plasmon/phonon  scattering,  band  filling,  quantum  mechanical  transmission,  and  band  bending 
are  studied  in  detail.  The  spectral  density  function  technique  is  also  discussed  in  conjunction 
with  the  possible  inclusion  of  collision  broadening,  a  quantum  mechanical  effect,  into  a  semiclas¬ 
sical  Monte  Carlo  simulation. 

The  nonlinear  transport  study  in  the  present  work  essentially  consists  of  two  aspects.  The 
first  topic  is  hot  electron  transport  in  GaAs,  focusing  on  the  electron  impact  ionization  effects. 
The  dependence  of  impact  ionization  rates  on  the  details  of  the  band  structure  is  investigated  by 
using  two  (local  and  nonlocal)  empirical  pseudopotential  methods.  The  spatial  evolution  of  the 
ionization  rate  and  the  average  electron  energy  are  studied  in  nonuniform  fields  characteristic  of 
p+-n  junctions  and  Schottky  barriers.  An  appropriate  definition  of  the  imapct  ionization  rate  in 
transient  problems  is  discussed  along  with  the  Keldysh  model.  The  effects  of  field  fluctuations 
due  to  the  random  distribution  of  dopants  are  studied  as  well,  and  are  found  to  have  only  lim¬ 
ited  importance  on  the  enhancement  of  the  impact  ionization  rate.  The  possibility  of  new  GaAs 
electron-emitting  diodes  is  explored  numerically  and  compared  with  the  corresponding  Si  devices 

The  second  aspect  deals  with  the  effects  of  conduction  band  discontinuities  on  the  electron 
transport.  In  particular,  one-dimensional  heterostructures  are  modeled  to  study  the  nonlinear 
transport  across  heterointerfaces.  First,  two  heterostructure  avalanche  photodiodes  are  studied 
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Particular  attention  is  paid  to  microscopic  features  such  as  the  evolution  of  the  energy  and 
momentum  distribution  function  in  space  and  time,  from  which  one  can  deduce  macroscopic 
results.  It  is  found  that  overheating,  enhanced  energy  relaxation,  and  carrier  confinement  as  a 
consequence  of  the  structure  in  real  space  have  a  pronouned  influence  on  the  energy  and  momen¬ 
tum  distribution.  As  a  result,  the  energy  distribution  can  have  a  structure  which  directly 
reveals  the  band  structure  (of  the  material).  The  dependence  of  the  impact  ionization  rate  on 
the  band  structures  of  neighboring  layers  is  also  addressed. 

The  effects  of  a  nonequilibrium  phonon  distribution  on  the  electron  transport  are  studied. 
The  phonon  distribution  can  be  perturbed  whenever  a  significant  number  of  phonons  (especially 
zone-center  phonons)  are  emitted  on  a  time  scale  shorter  than  the  phonon  decay  time  constant. 
Such  a  situation  can  occur  when  a  large  number  of  carriers  propagate  across  an  interface,  experi¬ 
ence  an  abrupt  energy  gain,  and  subsequently  relax  through  strong  phonon  emission.  An  algo¬ 
rithm  is  developed  for  the  microscopic  analysis  of  phonon  dynamics.  It  is  observed  that  the  hot 
phonons  change  the  scattering  rate  considerably  and  heat  the  electron  energy  distribution  due  to 
the  reabsorption  of  nonequilibrium  phonons. 

To  investigate  the  quasiballistic  nature  of  electron  motion  and  to  examine  the  validity  of 
hot  electron  spectroscopy,  the  tunneling  hot  electron  transfer  amplifier  structures  are  studied  at 
4.2  K.  Special  attention  is  paid  to  the  transport  in  the  collector  barrier,  a  region  in  which  the 
tunneling  hot  electron  transfer  amplifier  differs  from  the  similar  planar  doped  barrier  transistor 
Our  model  includes  the  effects  of  coupled  plasmon/phonon  interaction  and  the  Pauli  exclusion 
principle  along  with  more  conventional  features.  The  numerical  results  demonstrate  the 
existence  of  nearly  ballistic  transport  in  the  base  and  the  collector  barrier,  and  confirm  that  the 
experiments  can  indeed  measure  the  energy  distribution  of  injected  ballistic  electrons.  The  dev¬ 
ice  characteristics  such  as  transfer  ratio  and  transit  time  are  also  investigated  in  detail,  and  dis¬ 
cussed  along  with  the  optimum  operating  conditions. 
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